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ry^ We investigate the perturbative and nonperturbative renormalization of composite oper- 

ators in lattice QCD restricting ourselves to operators that are bilinear in the quark fields 
(quark- antiquark operators). These include operators which are relevant to the calcula- 
tion of moments of hadronic structure functions. The nonperturbative computations are 
based on Monte Carlo simulations with two flavors of clover fermions and utilize the Rome- 
Southampton method also known as the RI-MOM scheme. We compare the results of this 
approach with various estimates from lattice perturbation theory, in particular with recent 
Qh two-loop calculations. 
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The investigation of hadron structure has become a central topic of lattice QCD. In many 



CO 

cases this involves the evaluation of matrix elements of local operators between hadron states. For 
example, (moments of) generalized parton distributions can be extracted from matrix elements of 
quark-antiquark operators, i.e., operators composed of a quark field, its adjoint and a number of 
gluon fields entering through covariant derivatives which act on the quark fields. In general such 
operators have to be renormalized. In this process the operator of interest may receive contributions 



also from other operators, i.e., it may mix with these additional operators. On the lattice, mixing 
occurs more frequently than in the continuum due to the restricted space-time symmetries. Since 
in the end one wants to make contact with phenomenological studies, which almost exclusively 
refer to operators renormalized in the MS scheme of dimensional regularization, one needs the 
renormalization factors leading from the bare operators on the lattice to the MS operators in the 
continuum. 

The most straightforward approach towards the calculation of renormalization factors is based 
on lattice perturbation theory (for a review see Ref. [l|). Unfortunately, this method meets with 
some difficulties. First, perturbation theory on the lattice is computationally much more complex 
than in the continuum and therefore the calculations rarely extend beyond one- loop order (see, 
however, Refs. 0-0]). Second, lattice perturbation theory usually converges rather slowly so that 
the accuracy of perturbative renormalization factors is limited. Identifying one source of these 
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poor convergence properties, Lepage and Mackenzie proposed as a remedy the so-called tadpole 
improved perturbation theory [5[]. Still, considerable uncertainty remains. Third, mixing with 
operators of lower dimension cannot be treated by perturbation theory. 

In special cases, when the renormalization factors contain no ultraviolet divergences, a nonper- 
turbative determination is possible with the help of Ward identities 0]. However, there are many 
interesting operators that cannot be renormalized by this method. 

A general nonperturbative approach to renormalization has been developed within the 
Schrodinger functional scheme (see, e.g., Refs. 0, 0], reviews are given in Refs. BQ)- 

In this 

method the finite size of the lattices employed in the simulations (along with appropriate boundary 
conditions in Euclidean time) is used to set the renormalization scale. In the end continuum per- 
turbation theory is employed to convert the results from the Schrodinger functional scheme to the 
MS scheme. Though theoretically appealing the practical implementation of the procedure requires 
a lot of effort and has to be repeated for every new operator again from the very beginning. 

Another nonperturbative method for computing renormalization coefficients of arbitrary quark- 
antiquark operators is the Rome-Southampton method (also known as the RI-MOM scheme) in- 



troduced in Ref. 11]. It mimics the procedure used in continuum perturbation theory. The basic 



objects are quark two-point functions with an insertion of the operator under consideration at 
momentum zero. These are computed in a suitable gauge, e.g., the Landau gauge. In continuum 
perturbation theory the two-point functions are calculated order by order in an expansion in pow- 
ers of the strong coupling constant while in the Rome-Southampton method they are evaluated 
within a Monte Carlo simulation on the lattice. In order to extract from these data renormalization 
factors which yield renormalized operators in the MS scheme in the continuum limit one needs a 
renormalization condition which is applicable to lattice data as well as to perturbative continuum 



results. A suitable condition has been given in Ref. 11]. 

Compared with the Schrodinger functional approach the Rome-Southampton method is dis- 
tinguished by its relatively simple implementation. Furthermore, one can deal with all desired 
operators in a single simulation. On the other hand, the Schrodinger functional method is explic- 
itly gauge invariant, while the Rome-Southampton method requires gauge fixing. 

In a previous publication [l2| we have performed an extensive study of nonperturbative renor- 
malization for a variety of quark-antiquark operators using the Rome-Southampton method, mo- 
tivated by our investigations of hadron structure functions. This was done with Wilson fermions 
in the quenched approximation for two values of the lattice spacing. Later on, these studies have 
been extended to improved Wilson fermions, based on quenched simulations at three values of the 
lattice spacing [ljj. Meanwhile we are using gauge field configurations generated with two flavors 
of dynamical quarks, which has made a reconsideration of renormalization necessary. 

In this paper we present results for renormalization factors obtained within the Rome- 
Southampton approach with nj = 2 dynamical quarks. We work with nonperturbatively 0(a)- 
improved Wilson fermions (clover fermions). The operators are, however, not (yet) improved. 
We continue to apply the momentum sources introduced in Ref. [l3|. In addition, we have re- 
fined our approach, subtracting lattice artefacts through one-loop boosted perturbation theory. 
As in Ref. [12j we consider only flavor-nonsinglet quark-antiquark operators. For some thoughts 



concerning flavor-singlet operators see Refs. 14. 1 151] . 



The paper is organized as follows: After introducing in Sec. [XT] the operators to be studied 
we explain the method of nonperturbative renormalization in Sec. IHIi Our implementation of this 
method employing momentum sources is described in Sec. lIVI After a brief overview over our gauge 
field configurations in Sec. [V]we discuss the chiral extrapolation of our data in Sec. IVI1 Section l"VIII 
reviews formulae from continuum perturbation theory that will be needed in the analysis. Results 
from lattice perturbation theory are compiled in Sec. IVIIIl Section IIXI explains how we apply 
lattice perturbation theory in order to subtract lattice artefacts. In Sec. [X] we describe our method 
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of extracting the renormalization factors from the Monte Carlo data. The results (perturbative 
as well as nonperturbative) are then presented and discussed in Sec. IXII Finally, we present our 
conclusions in Sec. IXIII Some technical details are explained in the Appendices. 



II. THE OPERATORS 

In the Euclidean continuum we want to study the operators 



CW-**n = u^D^-'-D^d, (1) 
0^..., n = u W D tll ---D tMn d, (2) 
tfLm-Hn = ua^D^ ■ ■ ■ D^d (3) 



(with o^ v = (i/2)[7/i)7i/] an d = — D^) or rather 0(4) irreducible multiplets with definite 
charge conjugation parity. In particular, we obtain twist- 2 operators by symmetrizing the indices 
and subtracting the traces. We have given the quark fields definite flavors (assumed to be degen- 
erate) in order to make apparent that we are considering the flavor-nonsinglet case. Hence the 
twist-2 operators do not mix and are multiplicatively renormalizable. 

Working with Wilson fermions it is straightforward to write down lattice versions of the above 
operators. One simply replaces the continuum covariant derivative by its lattice analogue. How- 
ever, 0(4) being restricted to its finite subgroup H(4) (the hypercubic group) on the lattice, the 
constraints imposed by space-time symmetry are less stringent than in the continuum and the 



possibilities for mixing increase 16l4l9l] 



While the H(4) classification for operators Ouui —u n an d Cuui-u„ with n < 3 has been treated 
in detail in Ref. [18(], we have to refer to Ref. [201 ] for the operators OT^^.,. „ . Note however that 
the classification of the latter operators for n < 2 can be derived from the results presented in 



Ref. 



In our investigations of hadronic matrix elements we have considered the following operators 



whose renormalization factors have already been studied in Ref. [12 



(4) 
(5) 
(6) 
(7) 
(8) 
(9) 
(10) 

(11) 

Their labels refer to the structure function moments that they determine. These operators have 
been selected such that they have a definite transformation behavior under H(4) (i.e., belong to 
an irreducible multiplet) with as little mixing as possible. Moreover we have tried to minimize 
the number of nonzero momentum components required in the evaluation of the hadronic matrix 
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= 0{44} - 
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= ^{114} - 
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— CO 5 

— ^{124} ' 




Or 2>a 


= ^{14} ' 




Or 2 , b 


- U {44} 


3 (^{11} + ^{22} + ^{33}) ' 


O r3 


= ^{114} " 


" \ (^{224} + ^{334}) ' 



elements. Note, however, that in the numerical simulations reported in Ref. [13|] different operators 
(though from the same H(4) multiplets) have been used for the matrix elements V3 and V4. 

The operators V2a and V2b transform according to inequivalent representations of H(4), 
although they belong to the same irreducible 0(4) multiplet in the continuum. Therefore their 
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renormalization factors calculated on the lattice need not coincide. The same remark applies to 
O r2 a and O r2b . 

Since our matrix element calculations now involve additional operators, not considered in 
Ref. |12i | . we have extended the above list by the following operators, again guided by the H(4) 



classification given in Refs. [18|, [2 



(12) 
(13) 
(14) 
(15) 

+" ^{233} , (16) 
+" Cll{24} + ^42{14} ' ( 17 ) 
13{23} + ^23{13} ~~ ^41{24} ~ ®42{14} • ( 18 ) 

The operator V3a yields the same structure function moment as V3 . However, in contrast to 
V3 it cannot mix with any operator of the same or lower dimension. On the other hand, it has 
the disadvantage that one needs spatial momenta with two nonvanishing components in order to 
extract the moment v% from its matrix elements. The latter fact is the reason why we did not 
employ it in our previous investigations of nucleon structure. The operators constructed from OT. 
are relevant for transversity. 

Furthermore, we have studied the following operators without derivatives ("currents"): 



o V3 , a 


= C{124} > 




o hl , a 


- U l{23} ' 




o hl , b 


— u 122 


u 133 : 


o h2 , a 


— ^4{123} 




o h2 , b 


— U l{122} 


~ ^{133} " 


o h2 , c 


— ^13{23} 


+ ®23{13} ~ 


o h2 , d 


— ^1211 





O s = ud, (19) 
O p = u-ysd, (20) 



Ol = «7 M d, (21) 



= u W d, (22) 



A' 

0%, = ua^d, (23) 

where all quark fields are taken at the same lattice point. Finally we have also considered the 
quark wave function renormalization constant Z q . 

In Table [J we list all operators studied, along with the H(4) representation they belong to and 
their charge conjugation parity C. 

While in the evaluation of hadronic matrix elements the members of a given operator multiplet 
require different momentum components to be nonzero and hence are of different usefulness, such 
distinctions do not matter in our computation of renormalization factors. Therefore we consider 
not only individual operators but also complete operator bases for the representations under study. 
The representations studied and the chosen bases are given in Appendix lAl 

Concerning the mixing properties a few remarks are in order. Mixing with operators of equal 
or lower dimension is excluded for the operators V2 a , V2 b , V3 a , O a2 , O r2 a , O r - 2 b , a , Oh x b , 
Oh 2 a , Oh 2 b , Oh 2 c as well as for the currents. 

The case of the operator O v:i , for which there are two further operators with the same dimension 
and the same transformation behavior, is discussed in Refs. |l8l. Similarly, 0^ 2 d could mix 
with another operator of the same dimension. The operators Vi , O r3 , on the other hand, could 
in principle mix not only with operators of the same dimension but also with an operator of one 
dimension less constructed from O t ... A few more details on the mixing issue can be found in 
Ref. [ljj, in particular in Appendix B. 
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TABLE I: Operators and their transformation behavior [H|,[2(|. The charge conjugation parity is denoted 



by C. 
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Repr. 
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Op. 


Repr. 
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+1 


o hl , b 
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+ 1 


O p 


rP 


+1 


a 3 


r<»> 


-1 


ol 




-1 


o V3 , a 




-1 


o£ 




+1 


Or, 


T (8) 


+1 


o T » v 




-1 


O a2 


T (4) 


+1 


o V2 , a 




+1 


o h2 . a 


r (3) 


-1 


o V2 , 




+1 


o h2 , b 


T (3) 


-1 


O r2 , a 


r (6) 


-1 


o h2 , c 


r (6) 


-1 




r (3) 


-1 


Oh 2<d 




-1 


o hl , a 


T (8) 


+1 


o V4 




+1 



Our analysis ignores mixing completely. This seems to be justified for O vs . Here a perturbative 
calculation gives a rather small mixing coefficient for one of the mixing operators 17|, [19|], whereas 
the other candidate for mixing does not appear at all in a one-loop calculation of quark matrix 
elements at momentum transfer zero, because its Born term vanishes in forward direction. The 
same is true for all operators of dimension less or equal to 6 which transform identically to Vi : 
Their Born terms vanish in forward matrix elements, hence they do not show up in a one-loop 
calculation at vanishing momentum transfer. In the case of O r3 , however, the mixing with an 
operator of lower dimension is already visible at the one-loop level even in forward direction. 
Nevertheless, the nucleon matrix elements of the operators mixing with O m and O va seem to be 
small, at least in the quenched approximation 131 ]. 



III. THE METHOD 



We calculate our renormalization constants with the help of the procedure proposed by Mar- 
tinelli et al. [ll|] (the Rome-Southampton approach). It follows closely the definitions used in 
(continuum) perturbation theory. We work on a lattice of spacing a and volume V in Euclidean 
space. For a fixed gauge let 

G a p(p) = ^ E z-' lHx - y) Mx)0{z)d p {y)) (24) 

x,y,z 

denote the nonamputated quark-quark Green function with one insertion of the operator O at 
momentum zero. It is to be considered as a 12 x 12 matrix in the combined color and Dirac space. 
The corresponding vertex function (or amputated Green function) is given by 

T(p) = S- 1 (p)G(p)S- 1 (p), (25) 

where for q = u or q = d 

x,y 
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denotes the quark propagator. We define the renormalized vertex function by 

r R (p) = Z- x ZT{p) (27) 
and fix the renormalization constant Z by imposing the renormalization condition 

£tr (rR^rBorn^)" 1 ) = 1 (28) 

at p 2 = fj,p, where \i v is the renormalization scale. So Z can be calculated from the relation 

z^z^tx (r(p)rBorn(p)- 1 ) = 1 (29) 

with p 2 = /ip. Here TBom(p) is the Born term in the vertex function of O computed on the lattice, 
and Z q denotes the quark field renormalization constant. The latter is taken as 

7 f s _ tr H Ea 7a sm(ap x )aS- 1 (p)) 

again at p 2 = fj, 2 . Aiming at a mass- independent renormalization scheme we finally have to 
extrapolate the resulting values of Z to the chiral limit. 

Note that there will be no 0(a) lattice artefacts in Eq. (|29p. because they come with operators 
of opposite chirality in the vertex function, and these drop out when the trace is taken. Still, 
matrix elements of the renormalized operators will in general have 0(a) lattice artefacts because 
the operators are not improved. Once improved operators are available one can evaluate their 
renormalization factors using the methods described in this paper. 

Equations (j28[) and (|30p (in the chiral limit) together define a renormalization scheme of the 



momentum subtraction type which is called RI' — MOM scheme Here RI stands for "regular- 
ization independent" . This nomenclature refers to the fact that the definition of the RI' — MOM 
scheme does not depend on a particular regularization - here we have used a lattice cutoff just for 
definiteness and because the lattice regularization will be the basis of our numerical investigations. 
The MS scheme, on the other hand, can only be defined within dimensional regularization and is 
therefore restricted to perturbation theory. 

The RI' — MOM scheme differs from the RI-MOM scheme only in the definition of the quark 
field renormalization constant, which in the RI' — MOM scheme is more suitable for the numerical 
evaluation. 

In general, the RI' — MOM scheme will not agree with any of the momentum subtraction 
schemes used in continuum perturbation theory. It is therefore desirable to convert our results 
to a more popular scheme like the MS scheme. Another reason for converting to the MS scheme 
lies in the fact that many of the operators discussed in this paper appear in the operator prod- 
uct expansion along with the corresponding Wilson coefficients, which are generally given in the 
MS scheme. Hence we have to perform a finite renormalization leading us from the RI' — MOM 
scheme to the MS scheme if we want to use our renormalized operators together with the per- 
turbative Wilson coefficients. This finite renormalization factor can be computed in continuum 
perturbation theory using, e.g., dimensional regularization. The details needed for the evaluation 
of this factor will be discussed in Sec. I VIII 

If the operator under study belongs to an 0(4) multiplet of dimension greater than 1, i.e., 
if it carries at least one space-time index, the trace in Eq. (j29|) will in general depend on the 
direction of p. This has the immediate consequence that the renormalization condition (|28h violates 
0(4) covariance even in the continuum limit. In the continuum, one can restore 0(4) covariance 
by a summation over the members of the 0(4) multiplet. On the lattice, each operator when 
renormalized according to Eq. ([29p has in general its own Z factor, and only after conversion to a 
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covariant scheme all operators in an irreducible H(4) multiplet will have the same renormalization 
factor. However, it is also possible to define a common Z factor for all members of such an H(4) 
multiplet already in the RI-MOM framework by taking a suitable average. If j = 1,2, . . . , N labels 
the members of the chosen basis of the multiplet we can average over this basis and calculate Z 
from 

1 N 

z ^ z n £ ^ tr (rM^ipr 1 ) = i (3i) 

with p 2 = fj, 2 . This procedure has two advantages. It is simpler than working with a different 
renormalization factor for every single operator, and it leads to a smoother dependence of the 
results on p 2 , because it reduces the amount of 0(4) violation. 

The bases actually used in our calculations are given in Appendix Whenever we want to 
refer to this averaging procedure we shall write the respective operator with a bar on top, i.e., V2 a 
means precisely the operator @ while V2 a refers to a result for the operator multiplet (|Aip . and 
analogously for the currents. 

Ideally, the scale fi p at which our renormalization constants are defined should satisfy the 
conditions 

1/L 2 « A^ CD « fx 2 p « 1/a 2 (32) 

on a lattice with linear extent L. The inequality Aq CD <C fi 2 should ensure that we can safely 
use (continuum) perturbation theory to transform our results from one scheme to another. The 
inequality fi 2 <C 1/a 2 is supposed to keep discretization effects small. So we have to find a way 
between the Scylla of difficult to control nonperturbative effects and the Charybdis of lattice 
artefacts. Whether in a concrete calculation the conditions (|32p may be considered as fulfilled 
remains to be seen. 

Let us finally comment on our notation for the renormalization scale. In the case of a general 
scheme S we use the letter M, in the case of the MS scheme we use [i. We take fj, p when dealing 
with the Rl' — MOM scheme and fiM in the case of the MOM scheme to be defined below. 



IV. NUMERICAL IMPLEMENTATION 



Let us sketch the main ingredients of our calculational procedure 12]. To simplify the notation 
we set the lattice spacing a = 1 in this section. Moreover we suppress Dirac and color indices. 
In a first step the gauge field configurations are numerically fixed to some convenient gauge, the 
Landau gauge in our case [211 ] . Representing the operator under study in the form 

J20(z) = J2q(z)J(z,z')q(z') (33) 

z z,z' 

we calculate the nonamputated Green function (j24f) as the gauge field average of the quantity 

d ^ = V ^ z- lHx - y) S{x,z)J{z,z')S{z',y), (34) 

x,y,z,z' 

constructed from the quark propagator S on the same gauge field configuration. Working in the 
limit of exact isospin invariance we do not have to distinguish between u and d propagators. With 
the help of the relation 

S{x, y) = j 5 S(y , x) + 75 (35) 
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we rewrite G(p) as 



Gip) = 7 £ 75 (e s(z, x yA 15 j(z, z') s(z', y ypy) . 

z,z' \ x / \ y / 



(36) 



The quantities 



z,x)e ip - x (37) 



appearing in this expression can be calculated by solving the lattice Dirac equation with a momen- 
tum source: 

£ M(y, z) S(z, x^A = . ( 38 ) 

Here M(x, y) represents the fermion matrix. So the number of required matrix inversions is pro- 
portional to the number of momenta considered. But the quark propagators, which we need for 
the amputation and the computation of the quark wave function renormalization, are immediately 
obtained from the quantities already calculated. 

Strictly speaking, one should evaluate the quark propagators going into the calculation of S(p) 
in Eq. (126p on configurations that are statistically independent of those used for the computation 



of the Green functions ()24j) in order to avoid unwanted correlations 78J]. If we calculate these 
expectation values on two independent ensembles, statistical fluctuations in the quark propagators 
and the Green functions (I24D are uncorrelated, and (I3ip gives a good estimate of the true Z. If 
we calculate both expectation values on the same ensemble, the fluctuations will be correlated, 
particularly if the configuration number iV con f is small. It can be shown [22] that this will give a 
bias proportional to l/N con {. In our case statistical fluctuations are very small, because of our use 
of momentum sources, so we do not expect a large problem. Indeed, estimating the bias introduced 
by our procedure (see, e.g., Refs. [zj [23| for appropriate methods) we confirm this expectation. 
Note that such correlations exist also in the calculation of hadronic matrix elements from ratios of 
correlation functions. However, given the typical number of configurations used in these studies, 
the bias should be very small. 

Another computational strategy would be to choose a particular location for the operator. 
Translational invariance ensures that this will give the same expectation value after averaging 
over all gauge field configurations. For this method we need to solve the Dirac equation with a 
point source at the location of the operator and (in the case of extended operators) for a small 
number of point sources in the immediate neighbourhood. For operators with a small number of 
derivatives the point source method would require fewer inversions, but it turns out that relying 
on translational invariance increases the statistical errors. 

The required gauge fixing necessarily raises the question of the influence of Gribov copies. For- 
tunately, investigations of this problem indicate that the fluctuations induced by the Gribov copies 
are not overwhelmingly large_and may be less important than the ordinary statistical fluctua- 
tions 



2J,|2J] (see also Ref. [26j). 
Since the numerical effort is proportional to the number of momenta, the proper choice of the 
momenta considered is of particular importance. In order to minimize cut-off effects we choose 



them close to the diagonal of the Brillouin zone and achieve for most operators an essentially 
smooth dependence on the renormalization scale fj,p. It goes without saying that in this way we 
cannot eliminate lattice artefacts completely. However, more sophisticated strategies for coping 



with the cut-off effects such as those suggested in Refs. [27], l28|] would require the use of far more 
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TABLE II: Simulation parameters j3, n = K sca , clover coefficient Csw and lattice volume along with the 
corresponding values of the pion mass in lattice units. 



p 


K 


C SW 


V 


am, 


5.20 


0.1342 


2.0171 


16 3 x 32 


0.5847(12) 


5.20 


0.1350 


2.0171 


16 3 x 32 


0.4148(13) 


5.20 


0.1355 


2.0171 


16 3 x 32 


0.2907(15) 


5.25 


0.1346 


1.9603 


16 3 x 32 


0.4932(10) 


5.25 


0.1352 


1.9603 


16 3 x 32 


0.3821(13) 


5.25 


0.13575 


1.9603 


24 3 x 48 


0.25556(55) 


5.25 


0.1360 


1.9603 


24 3 x 48 


0.18396(56) 


5.29 


0.1340 


1.9192 


16 3 x 32 


0.5767(11) 


5.29 


0.1350 


1.9192 


16 3 x 32 


0.42057(92) 


5.29 


0.1355 


1.9192 


24 3 x 48 


0.32696(64) 


5.29 


0.1359 


1.9192 


24 3 x 48 


0.23997(47) 


5.29 


0.1362 


1.9192 


24 3 x 48 


0.15784(75) 


5.40 


0.1350 


1.8228 


24 3 x 48 


0.40301(43) 


5.40 


0.1356 


1.8228 


24 3 x 48 


0.31232(67) 


5.40 


0.1361 


1.8228 


24 3 x 48 


0.22081(72) 


5.40 


0.13625 


1.8228 


24 3 x 48 


0.19053(62) 


5.40 


0.1364 


1.8228 


24 3 x 48 


0.1535(13) 



momenta than we can afford when working with momentum sources. As the treatment of lattice 
artefacts is a subtle issue anyway we have decided to keep the advantage of small statistical errors 
provided by the above procedure and to deal with the discretization errors in a different manner. 

A specific lattice artefact is caused by the 0(a) chiral symmetry breaking term of the quark 
propagator. In position space this term is concentrated at very short distances, for fermion actions 
obeying the Ginsparg- Wilson condition it is even exactly a delta function. In momentum space 
it gives rise to the Wilson mass term, ~ ap 2 in the inverse propagator. The authors of Ref. [2^] 
discuss one method of suppressing this artefact (for an earlier treatment of the same effect see 
Ref. [30]). Here, our approach to this problem is to define Z q from Eq. (|30p . in which the trace 
removes the Wilson mass term, and to suppress the remaining 0(a 2 ) lattice artefacts by using the 
perturbative subtraction scheme described in Sec. IIXI 



V. MONTE CARLO ENSEMBLES 



We work with two degenerate flavors of nonperturbatively improved Wilson fermions (clover 
fermions). For the explicit form of the fermionic action see, e.g., Ref. 30J] . As our gauge field action 
we take Wilson's plaquette action. In Table HI] we collect the parameters of our simulations, f3, 
K = K seau , the clover coefficient c sw and the lattice volume along with am^, the pion mass in lattice 
units. Table IIIII contains the critical hopping parameters k c as well as the values of the Sommer 
parameter r^/a and the average plaquette P = (htrUn) in the chiral limit [3l[ which are employed 
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TABLE III: Critical hopping parameters k c along with chirally extrapolated values for the Sommer parameter 
r 1 a and the average plaquette P. 



(3 n c r /a P 

5.20 0.136008(15) 5.454(59) 0.538608(49) 

5.25 0.136250(7) 5.880(26) 0.544780(89) 

5.29 0.136410(9) 6.201(25) 0.549877(109) 

5.40 0.136690(22) 6.946(44) 0.562499(46) 




FIG. 1: Z m '~ MOM for the operator O t at /3 = 5.29, k = 0.1362 on a 24 3 x 48 lattice. The curves represent 
splines with two interior knots fitted to the data and to the data ± the statistical error. 



in this paper. Note that the results for the chiral extrapolation of ro/a given here are based on a 



larger set of data than that used in Ref. 31]. 



The statistical errors will be calculated by means of the jackknife procedure. 



VI. CHIRAL EXTRAPOLATION 

As already mentioned in Sec. IIIII we have to extrapolate our results obtained at nonvanishing 
quark masses to the chiral limit. This will be done for each (3 at fixed values of p 2 . In the cases 
where the simulations for different values of k have been performed on different volumes, i.e., for 
(3 = 5.25 and f3 = 5.29, the sets of momenta used depend on k, and some kind of interpolation is 
required. For this purpose we fit the data on the larger lattices (24 3 x 48) with cubic splines in 
ln(a?p 2 ). Except for very small momenta, which will not influence the final results, these fits yield 
a very good description of our data. An example is shown in Fig. [TJ 

Of course, "wiggles" in the data (caused by lattice artefacts) will be smoothed out by this 
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interpolation. These wiggles are less pronounced on the larger lattices than on the smaller ones. 
That is why we have chosen to work with the momenta coming from the smaller lattices so that 
we can use the data on these lattices directly without any interpolation and have to interpolate 
only the results obtained on the larger lattices. 

Alternatively, one could use the interpolation for all k values. This however leads to negligible 
differences in the final results. 

For the chosen momenta we can then extrapolate our data to the chiral limit. This is done 
linearly in (rom n ) 2 , i.e., by means of a fit of the form 

Z = z + z 1 (r m n ) 2 , (39) 

where the fit parameter z$ is identified with the desired value of the renormalization factor in the 
chiral limit. Note that this is essentially a linear fit in the quark mass. The ansatz is motivated 
by the fact that in perturbation theory the leading quark mass dependence is linear as the chiral 
limit is approached (see, e.g., Ref. [30]). 

With the possible exception of the smallest momenta these fits work well. Examples are shown 
in Figs. [21 [3l Nevertheless, we have also performed quadratic extrapolations of the form 

Z = z + z 1 (r m 7T ) 2 + z 2 (r Q m 7T ) 4 ' (40) 

in order to get an idea of the impact of the chiral extrapolation on the final results. (Note that 
this corresponds to a three-parameter fit at three data points for /3 = 5.20.) 
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FIG. 2: Chiral extrapolation for O a2 at j3 = 5.40. 

However, there is an exceptional case where these simple extrapolations are not trustworthy. 
This is the pseudoscalar density O p . In this case one expects that Z = Zp vanishes with the quark 
mass rriq, because Zp 1 develops a pole in m q . Therefore we follow Ref. [32| and try to subtract the 
pole contribution using a fit of the form 

— = s — \-si + s 2 am q (41) 

Zp am q 
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FIG. 3: Chiral extrapolation for O s (subtracted data, as explained in Sec. IIXI) at (3 — 5.40. 



with 



- - 

2\k 



(42) 



Here k c is the critical hopping parameter defined for fixed /3 by the vanishing of the pseudoscalar 
mass m^. The values used in this paper can be found in Table [TTT1 The fit parameter s\ is then 
identified with the inverse of Zp in the chiral limit. Examples of such fits are shown in Fig. [U The 
curvature in the data is clearly visible establishing the existence of the Goldstone pole. 

How can we judge the reliability of the resulting numbers? From the operator product expan- 
sion (33j,|34j we expect that so is inversely proportional to /jz, i.e., /jz so should become independent 



of Hp. Therefore we plot a 2 fip so versus up 1 in Fig. [5j The fi 2 independence seems to be satisfied 
with reasonable accuracy. Thus we are confident that our extrapolation for Zp works fairly well. 
Nevertheless, the results for Zp must be considered with some caution. 



VII. INPUT FROM CONTINUUM PERTURBATION THEORY 



In Sec. IIIII we have explained how one can compute nonperturbative renormalization factors 
leading us from the bare lattice operators (at lattice spacing a) to renormalized operators in the 
RI' — MOM scheme (renormalization scale fi p ). In this section we collect results from continuum 
perturbation theory which will be needed for the conversion to standard renormalization schemes 
such as the MS scheme. 

If the operator O under study is multiplicatively renormalizable the operator renormalized in 
some scheme S at the scale M is related to the bare lattice operator Obare by 



O s (M) = Z^ arc (M,a)O bare 



(43) 
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FIG. 4: Chiral extrapolation of 1/Zp (subtracted data, as explained in Sec. IIX|) at j3 = 5.40. The symbols 
at am = represent the chirally extrapolated values, i.e., the quantity s±. 
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FIG. 5: a 2 fip so (from a fit to subtracted data, as explained in Sec. IDC [) as a function of (i* for the pseudoscalar 
density at /3 = 5.40. 

The scale dependence of the renormalized operator is determined by the anomalous dimension 



7' 



s = -M—lnZt 



dM 



bare ■ 



(44) 



Here the derivative is to be taken at fixed bare parameters, and it is implicitly assumed that the 
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cutoff has been removed in the end. In perturbation theory j s is expanded in powers of some 
renormalized coupling constant g s (M): 

7 =70 1g , +7i 1g , +7 2 1fl o +7 3 1fi o +••• (45) 



16vr 2 V !6vr 2 7 2 V ^vr 2 J 6 \ 16vr 2 

Note that the one- loop coefficient 70 is scheme independent. 

Similarly we define the quark field renormalization constant Z^ bare (M, a) in the scheme S so 
that the renormalized quark propagator is given by Zf h3Xe (M, a)5 b are- In the RI' - MOM scheme, 

~^ IOM is then specified by Eq. ([30]) or its continuum analogue. For the anomalous dimension 
of the quark field we adopt the definition 

7? = -M±\nZf ihaze . (46) 

The running of the coupling constant g (M) as the scale M is varied is controlled by the j3 
function 

P S = M^Lg s (M). (47) 

Again, the derivative is to be taken at fixed bare parameters and it is implicitly assumed that the 
cutoff has been removed in the end. The perturbative expansion of the (3 function can be written 

as 

s g s {Mf g s (Mf 5 5_W og gfCjg , (m 

In this case the first two coefficients (3q and (3\ are scheme independent. 
Integrating Eq. (|47jl we obtain 

M (J^„s (M ^^ ( 1 167f2 



A 5 "" Vl6vr 2y K ' ) 1 \2f3 g s {M) 2 



-xpW tf[*+±.™"L)\ (49) 



with the A parameter A s appearing as an integration constant. 

In the same spirit we define the so-called RGI (renormalization group invariant) operator, which 
is independent of scale and scheme, by 



Q RGi = AZ 5 (M)C 5 (M) = Z RGI (a)O b are (50) 



with 



and 



Z RGI (a) = AZ s (M)Zt re (M, «) > (52) 



where Z RGI depends only on a (or on the bare coupling parameter f3). Once we know RGI (or 
equivalently Z RGI ), multiplication with AZ S (M)~ 1 will allow us to evaluate Z^ re (M , a) and hence 
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the operator O (or rather its matrix elements) in any scheme and at any scale we like, provided 
we know the /3 and 7 functions sufficiently well. 

In the two-loop approximation, i.e., setting = 7^ = for n > 2, one can evaluate the integral 
in Eq. (fSTj) easily: 

Whenever we need them we evaluate the integrals in Eqs. (|49p and ()51[) exactly (by numerical 
methods) and do not re-expand them in g s (M). 

With the help of the methods described in Sees. [Ill] and [IV] we can compute Z^~ M (/j, p ) 
numerically for some range of scales fi p . Knowledge of AZ RI ~ MOM (n p ) would then permit us to 
compute Z RGI . However, being not covariant for most operators, the RI' — MOM scheme is not very 
suitable for evaluating anomalous dimensions. Therefore we will adopt a two-step procedure for 
computing Z RGl . In the first step we transform the numerical results for Z^ r ~ MOM to a covariant 
"intermediate" scheme S, e.g., the MS scheme, and in the second step we use the anomalous 
dimension and the /3 function in this scheme to compute AZ S and hence Z RGI . Thus we could in 
principle compute Z RGI (a) as 

Z RGI (a) = AZ S (M = v P )Zi^ MOM (M = H P )Z^ MOM (^ p , a) , (54) 

where ^/_ M q M denotes the finite renormalization factor leading from the RI' — MOM scheme to 
the scheme S and all the scales have been identified with the scale \i v initially set in the RI' — MOM 
scheme. 

The most obvious choice for S is of course the MS scheme. However, it will turn out to 
be advantageous to consider also a kind of (perturbative) momentum subtraction scheme, which 
we call MOM scheme. This is defined by requiring that in the renormalized vertex function the 
coefficient of the tree-level (or Born) term equals one at the renormalization scale \im- To make this 
definition unambiguous one has to specify in each case the basis used for the other contributions. 
The quark field renormalization constant is taken to be the same as in the RI' - MOM scheme: 

7MOM/ „\ _ 7 RI'-MOM/,, \ 

Z g,bare(^M,aJ - Z 9 bare [fiM,a). (55) 

The MOM scheme is covariant and rather "close" to the RI' — MOM scheme so that the conversion 
factor Z^\ ou from RI' - MOM to MOM usually differs less from one than the factor ^L M OM 
leading from RI' - MOM to MS. 

We can expand the conversion factor leading from RI' — MOM to MS in powers of the MS 
coupling constant: 



Z RI'-MOM(^) - 1 + C l 167r 2 + C 2 I 167r 2 I + C 3 I 16 ^ 2 I r ' ' ' ( 56 ) 



.MS 



Using 



z ™ (»)-i + c > 9MS ^ 2 + c' h MS ^ 2 ) , c > h MS M 2 ) .... (57) 
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we obtain for the conversion factor from RI' — MOM to MOM: 
7 MOM / \ _ f^RJ^MOMW 

^ri'-momW - m 

MOM^"' 

+ (Cl " Cl) 167T 2 + ( ° 2 " ClCl + Cl " ° 2) 1 167T 2 J 

+ (c 3 - C2C ; + C1C ; 2 - C1C ' 2 - c ; 3 + 244 - 4) + 



(58) 

If the vertex function r(p) is proportional to the Born term rBom(p) (as it happens in simple 
cases), the RI' — MOM scheme and the MOM scheme do not differ. So we have 

-^Rr-MOM^) = ^mom(m) (59) 

and consequently Cj = c-. However, in the generic case the matrix T(p) will contain also contribu- 
tions that are not a multiple of rBorn(p)- As we consider only operators which are multiplicatively 
renormalizable in the continuum these additional contributions are finite, but they make -^rj?_mom 
different from Z^q M and are responsible for the dependence of ^?_mom on ^ ne direction of the 
momentum p. An explicit example is discussed in Appendix [Bl 

Working with the MOM scheme it is quite natural to expand in a coupling constant which is 
similarly defined through a momentum subtraction procedure. Therefore we have also considered 



the MOMgg scheme as defined in Ref . 35|] . The corresponding coupling constant g MOM ss [ s related 
to the MS coupling constant g by 

16^ =^6^ + dl [-l6^) +d2 \-m^) +d3 [-w^) + '"' (60) 

where in the Landau gauge 
70 22 

516217 153 , /8125 4 \ 376 2 
d2 = ~m~ " - C3 " {-W + 3 C3 ) ^ + ' (62) 

304676635 299961 t 81825 . / 13203725 13339 . 1885 



C3 7T-C5 + 1onR + ~^^3 + -S-C5 "/ 



+ 



6912 64 ,J 64 w V 1296 27 

/ 580495 40 . \ 2 5680 



+ — C3 nj ni . (63) 

V 972 9 Sd 7 f 729 f V 7 



As usual, C„ denotes the value of Riemann's £ function at the argument n. Choosing S = MOM 
in Eq. (|51|) we shall always work with the MOMgg coupling. 

Using the above expressions, all expansions in the MS coupling may be rewritten as expansions in 
powers of the MOMgg coupling. For example, the coefficients ~f} AOM of the anomalous dimension 
in the MOM scheme expanded in powers of g MOM &s are related to the coefficients 7^ of the 
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anomalous dimension in the MS scheme expanded in powers of g MS by 



^MOM 

MOM 
72 

„,MOM 
73 



7l MS - 2A)Ci " ^170 , 



A1S 



2/Vi - 2/3 (24 " cf ) - 2di ( 7l MS - 2ci/3 ) + 7o (2d? - d 2 ) , 



72" 

7p-2^ vlo c' 1 ^;! ( / l 



(64) 
(65) 



MS j 



72 



„'2\ 



2/3i (2c' 2 - c' 



'2\ 



2 C ;/3i - 2/3 (24 - Cl 
2/3 (3c' 3 - 3c[c' 2 + e?) + 70 (5did 2 



.MS 



+ -2^4 (5d?-2d 2 ) 



5d? 



d 3 ) • 



(66) 



For the actual evaluation of the expansion coefficients q and c[ one starts from the bare vertex 
function for the operator under consideration computed in dimensional regularization and imposes 
the respective renormalization conditions yielding the renormalized vertex functions in the different 
schemes. These differ only by (scale dependent) factors, from which the desired conversion factors 
can be derived immediately, once the quark field renormalization factor has been extracted from 
the quark propagator. The results for the coefficients c±, c 2 , C3 and well as the sources 

from which we have taken the required perturbative vertex functions are given in Appendix O 



VIII. LATTICE PERTURBATION THEORY 

As long as no mixing with operators of lower dimension is involved it is possible to compute 
renormalization factors in lattice perturbation theory. Although straightforward in principle, the 
actual calculations tend to become rather cumbersome in practice. Hence they rarely extend beyond 
one-loop order (see, however, Refs. @~Sj])- This is a severe limitation since lattice perturbation 
theory converges rather slowly in most cases of interest. Therefore various improvement schemes 
have been devised, such as boosted perturbation theory and tadpole improvement jH]. 

In spite of these problems we want to compare our nonperturbative results with the corre- 
sponding values obtained in (improved) lattice perturbation theory. For the renormalization factor 
■^bareC^i a ) a straightforward application of one-loop lattice perturbation theory yields results of 
the form 

ZJ2,(M, a) pert = 1-^(70 ln(a M ) + C F A) , (67) 

where A = A(c S w) is a finite constant depending on the details of the lattice action and we have 
Cf = 4/3 for the gauge group SU(3). If mixing occurs, the single renormalization factor of a 
multiplicatively renormalizable operator is replaced by a matrix of Z factors. However, we shall 
neglect this complication and restrict ourselves to the matrix element on the diagonal corresponding 
to the operator under consideration. Working with an anticommuting 75 also in the continuum 
part of the calculation one obtains the values given in Table IIVI They do not depend on the 
particular operator but only on the H(4) multiplet to which the operator belon gs. In the case 
c sw = results for V4 , Oh lb , Oh 2b an d have already been obtained in Ref. [36J. Note that 
A is gauge invariant for our quark-antiquark operators, however, in the case of the quark wave 
function renormalization constant Z q it is not. The result given in Table ITVl corresponds to the 
Landau gauge. 

In order to obtain the corresponding results in tadpole improved perturbation theory we write 
(with /i = 1/a) for an operator with no covariant derivatives 
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TABLE IV: Finite contributions to the renormalization factors in lattice perturbation theory. The 

result for the quark field is given in the Landau gauge. 



Op. 


A(c sw ) 


Ref. 


s 


12.95241 + 7.73792c sw - 1.38038c| w 


[30] 


p 


22.59544 - 2.24887c sw + 2.03602c| w 


[30] 


ol 


20.61780 - 4.74556c sw - 0.54317c| w 


[30] 




15.79628 + 0.24783c sw - 2.25137c§ w 


[30] 




17.01808 - 3.91333c sw - 1.97230c| w 


[30] 




16.64441 - 2.24887c sw - 1.39727c| w 


[30] 


o V2>a 


1.27959 - 3.87297c sw - 0.67826cf w 


[30] 


o V2 . b 


2.56184 - 3.96980c sw - 1.03973cg w 


[30] 


Or 2 , a 


0.34512 - 1.35931c sw - 1.89255c§ w 


[30] 




0.16738 - 1.24953c sw - 1.99804c§ w 


[30] 


o hl , a 


1.25245 - 3.10180c sw - 1.59023c| w 


[37] 


o hl , b 


0.52246 - 2.99849c sw - 1.46224c| w 


[37J 


a 3 


-12.12740 - 2.92169c sw - 0.98166<| w 


[13,37,38] 


o V3 , a 


-11.56318 - 2.89800c sw - 0.98387c| w 


[37] 


Or, 


-12.86094 - 1.49316c sw - 1.68673c^ w 


[37J 


Oa 2 


-12.11715 - 1.51925c sw - 1.71846c| w 


[37J 


Oh 2 , a 


-11.54826 - 2.41077c sw - 1.51175c§ w 


[37] 


o h2 , b 


-11.86877 - 2.30651c sw - 1.34908cl w 


[37J 


Oh 2 , c 


-11.74773 - 2.36201c sw - 1.45084c| w 


[37] 


o h2 , d 


-12.9268 - 2.38849c sw - 1.3900c| w 


this work 




-25.50303 - 2.41788c sw - 1.12826c| w 


[13, 38] 


where 
and 


n = (itr%)3 = 1 - r£^C F 7r 2 + 0(g A ) 
A = A + (n D - l)vr 2 . 


(69) 
(70) 



This reflects the fact that one has operator tadpole diagrams and one leg tadpole diagram, which 
are of the same magnitude but contribute with opposite sign. It remains to make a physically 
reasonable choice for the expansion parameter ^lat- Here we identify <?lat with the boosted 
coupling 



(71) 
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Now we have two options. Either we stay with the expression (|67p and its tadpole improved 
analogue 



„2 

7MS / \ . _ 1-ric 



l-^hoHav) + C F A) 



(72) 



or we apply these formulae only at a fixed scale /i = (e.g., fio = 1/a) using the renormalization 
group to change fi: 

2£™(M, «)pe?t = A^)- 1 AZ M ( W )C( W , a) pert , (73) 
Z^M< G = A^^AZ^olC^^),. (74) 

The latter option seems preferable leading to the estimates AZ MS (/i )Z^ c (/io, a) pe rt and 

AZ MS (/Uo)^bare(/ i Oi °)ti for Z RGI . Working in the chiral limit we compute uo from the chirally 
extrapolated values for P = Mq given in Table Hill To be consistent with lowest order perturbation 
theory we set c sw = 1. 

Further improvement can be attempted by tadpole-improved, renormalization-group-improved 



boosted perturbation theory or TRB perturbation theory lj, |39j. This works as follows. In 
Eq. (144)) we have defined the anomalous dimension 7^ by differentiating Z^ axe with respect to 
the renormalization scale M at fixed cutoff and bare parameters. Alternatively one can keep the 
renormalized quantities fixed and take the derivative with respect to the cutoff, the lattice spacing 
o in our case. Then one obtains 

7 LAT = -aAi nZb 5 are . (75) 

Note that the derivative with respect to a also acts on c sw , unless only the tree-level value c sw = 1 
is used. The anomalous dimension 

7 LAT is 

to be considered ELS £1 function of some bare coupling 
constant #lat = Slat (a)- This could be the usual bare coupling g, but for our purposes it will be 
more advantageous to work with the boosted coupling <?□. Expanding in 5lat an d recalling that 
the one-loop coefficient 70 is universal we can write 

7" T («.AT)=7of|^ + 7r(f|f) 2 + - (76) 

Similarly we define 

oLATV \_ d^LAT o AlAT o ffLAT , n r Jl \ / 77 \ 

(3 {gLAT )_- a —-_-p —-p 1 ^—^ + O(g LAT ), (77) 

where (3q and j3\ have the same values as in the /3 function (|47p . Expressing Z RGI (a) = 
AZ s (M)Z^ are (M, a) in terms of 7 LAT and /3 LAT we find 



70 



,f|ff)' M %xnjr AT d„„r^ LA >> + ^ 

In the two-loop approximation we obtain 

vRGI _ ( ? o ti AT A , Pi A . . 

z -\ 2Po m?) y l + J Q m?) ■ (79) 



-YQ01 -7^ AT /3q 
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Choosing <7lat = 9a one has 

7i LAT = 7i° = 7P + 2p C F A(c sw ) + 167T 2 7 o (*i - \C F ) , (80) 



where [4CH43I] 

t x = 0.4682013 - (0.0066960 - 0.0050467c sw + 0.0298435c| w ) n f . (81) 
Tadpole improvement finally yields the result in TRB perturbation theory: 

70 „ 7001-71 fip | 2 i-i „ \ri Po 

Z KG1 _ l-n D ( 2B 9l Y*K( 1+ ft 9a \ {1 - nD)C ^ 

^o^j [ 1 + J T6^) ■ (82) 

Applying this formula we shall set again c sw = 1 to be consistent with lowest order perturbation 
theory. 

For the operators O s , O p , O^, O^, and 0^ u without derivatives two-loop calculations of the 
renormalization factors in lattice perturbation theory have recently appeared @>0|- The var i° us 
improvement schemes can be applied also to these two-loop expressions. However, the resulting 
formulae become considerably more complicated. Therefore we defer the corresponding discussion 
to Appendix |Pl 



IX. PERTURB ATIVE SUBTRACTION OF LATTICE ARTEFACTS 

In the perturbative form (|67|) of the renormalization factors the lattice spacing a only appears 
in the logarithm (and implicitly in the bare gauge coupling g). In the remaining contributions the 
limit a — > has been performed at fixed fi leading to the finite constant A. In this way all lattice 
artefacts vanishing like powers of a have been eliminated. However, there is no need to do so. In 
fact, a[i is not necessarily small in our Monte Carlo results. Hence it is worthwhile to keep a finite 
and to compare the lattice artefacts in the perturbative expressions with their nonperturbative 
counterparts. 

To do this we simply write down the one-loop integrals for general external momentum p and 
perform the integrations numerically. The integrals can no longer be reduced to a small number of 
standard integrals, they have to be done independently at each value of p and can only be obtained 
in numerical form. They will in general not only depend on p 2 but also on the direction of the 
momentum. 

For general p we can write the one-loop expression for Z?^ T ~ M M in the form 

^ e - M ° M (P> a) = 1 + r£^C F F(p, a) + 0(g 4 ) , (83) 

where the quark mass has been set equal to zero. Neglecting all contributions which vanish as 
a — > we get from F(p,a) the expression F(p,a); e.g., for the scalar density O s it is given by 

F(p,a) = 31n(aV) - 16.9524 - 7.73792c sw + 1.38038c| w (84) 

in the Landau gauge. The difference between F and F represents the lattice artefacts in one-loop 
perturbation theory. Though being 0(a 2 ), F — F can be fairly large for the momenta in the actual 
simulations. An example for the case of the scalar density O s is shown in Fig. [6l 

We can (and will) use this calculated difference to correct for the discretization errors in our 
lattice data (see also Ref. [3|). We take 



D(p,a) = F(p,a) - F(p,a) 



(85) 
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F(p.a) 




10 2 10 1 1 10 1 2 

a p 



FIG. 6: Lattice artefacts for the scalar density. The dotted straight line shows F(p, a), while the two other 
curves represent F(p, a) for the momentum directions indicated in the plot. The black squares denote the 
values of F(p, a) at the momenta used on the 16 3 x 32 lattices, which have been chosen close to the diagonal 
of the Brillouin zone. 

as an estimate of the perturbative discretization errors in our Monte Carlo renormalization con- 
stants Z h ~ M (p,cl)mc and define subtracted renormalization constants by 

^RI'-MOM/ \ ^RI'-MOM/ \ 9n n n / \ / CR \ 

Z bare IP, ajMC,sub = ^ bare (p, CtjMC - Jg^2 C F D {P, «) (86) 

employing boosted perturbation theory with u$ in gu taken at the respective k value, i.e., before the 
chiral extrapolation. Working consistently with one-loop perturbation theory we set c sw = 1. This 
procedure removes all the g 2 a 2 discretization terms in Z^~ M M (p, o)mc> leaving lattice artefacts 
0(g 4 a 2 ). As we shall see, the use of the boosted coupling g n seems to do a reasonable job of 
estimating the higher-order discretization effects. 

Unfortunately, this procedure becomes rather cumbersome for operators with more than one 
derivative. So we can use it only for the quark wave function renormalization, the currents and 
the operators with one derivative. 

X. EXTRACTING THE RENORMALIZATION FACTORS 

The simplest procedure for obtaining a value of Z RGI (a) would be to plot the right-hand side 
of Eq. ((SU), i.e., of the relation 

Z RGI (a) = AZ S (M = ^)Z| r _ MOM (M = ^4;™^, a) , (87) 

versus \x v and to read off Z RGI (a) in an interval of \x v where the inequalities (I32p are satisfied. In 
this region the value of Z RGI (a) would be independent of /j, p , i.e., one would observe a plateau, 
and one could determine the final result by fitting a constant to the data for Z RGI (a). Examples 
of such plots before and after the perturbative subtraction of lattice artefacts are shown in Fig. [71 
and in Fig. [8] subtracted and unsubtracted results for the quark field renormalization constant Z q 
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FIG. 7: Z RGI for the operator Ot before (upper plot) and after (lower plot) the perturbative subtraction 
of lattice artefacts. 

are directly compared at /3 = 5.20 and (3 = 5.40. Equivalently one could fit the values obtained for 
Z Le^ a ) = Z m>-MOM^p)Z^- MOM {n P ,a) in the plateau region by AZ S {^y 1 Z KG \a) (with 
Z RG1 (a) as fit parameter) using the most accurate perturbative expressions for ^rt'-mom an< ^ 
AZ S . 

However, in our actual simulations it is not so clear how well the inequalities (|32[) are fulfilled, 
and there are two effects to be considered that jeopardize the reliability of this approach. Firstly, 



there will be lattice artefacts which vanish like powers (up to logarithms) of a for a — > 271. |28I| . In 



order to reduce the corresponding contamination one would like to perform the fit at small values 
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FIG. 8: Z RGI for the quark field renormalization constant Z q before and after the perturbative subtraction 
of lattice artefacts. 

of Hp. Secondly, the truncation of the perturbative expansions in AZ S and ^rt'-mom wm produce 
noticeable effects in the region of small fi p leading, in particular, to a dependence of the results on 
the intermediate scheme S. In order to minimize the related uncertainties one would like to move 
the fit interval to large values of Hp- Because of these conflicting requirements it is a nontrivial 
matter to extract a final value for Z RGI from the data. 

Let us first investigate to which extent we can separate truncation effects from lattice artefacts. 
According to Eq. (JMJ), Z^ ie {fJ- P , a) = Z^ v _ MOM (fj, p )Z^ T ~ MOM (fi p ,a) can be written as 

Z^ tc (h p , a) = AZ 5 (^)- 1 Z RGI (a) . (88) 

In this way, we have factorized the dependence of Z^ ave (H p , a) on the renormalization scale Hp and 
on the cutoff a. Consequently we can write 

"^bare(^p; a ) = ^RGI \j\ ^bare 0^ a ') ■ ( 89 ) 

Hence multiplication by an appropriate (fjL p independent) scaling factor should bring the values 
of ^ifare(^P' a ) obtained for different values of a (or /3) onto a single curve representing a function 
f s (fi p ) of Hp only, provided \ip is small enough so that lattice artefacts can be neglected. Note that 
the ratio Z^ am (s/j. p , a) / 'Z^ &re (/j, p , a) for some fixed value of s, the so-called step scaling function, has 
a decent continuum limit. This coincides with f S {sHp) / f S {Hp) in the region where f s {Hp) is wen 
defined. 

In most cases, this collapse onto a single function works quite well for a reasonable range of 
renormalization scales, even if mixing is allowed. For an example see Fig. [9l So the factorization 
of Hp dependence and a dependence seems to be possible (except for the highest values of Hp)- 

However, the available perturbative results cannot describe the Hp dependence below Hp ~ 
5 GeV 2 , as exemplified by Fig. [7J It would be interesting to investigate whether this observation 
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FIG. 9: ^baro (perturbatively subtracted) for the operator Ot as a function of the renormalization scale 
fi p . The upper plot shows the actual results, in the lower plot they have been multiplied by suitable (/Lt p 
independent) scaling factors. 



is related to the claim 45] that DGLAP evolution cannot be used below Q 2 = 5 GeV 2 (at least 
in the region of larger values of Bjorken's variable x [7^]). Note also that an even later onset of 
(three-loop) perturbative behavior has been found for the quenched gluon propag ator 0,53- 

In Fig. [7]we have taken S = MOM and we have exploited the freedom to select the scheme S' for 
the coupling used in the perturbative expansion of %(_ MO m choosing S' = MOMgg. Generally, 
the plateaus in Z RGI look better for S = MOM than for S = MS. This may be due to the fact 
that the perturbative expansion of ^_ MO m seems to be better behaved for S = MOM. For S', 
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FIG. 10: Z RGl (perturbatively subtracted) for the operator 0\ lx a at /5 = 5.40 as a function of the renor- 
malization scale \x v . The open circles (filled squares) have been obtained with S = S' = MS (S — MOM, 
S' = MOMgg). 

the choice S' = MOMgg turns out to be preferable. An example comparing results obtained with 
5 = MS, S' = MS to results obtained with S = MOM, S' = MOMgg is shown in Fig. [KH Note 
that the difference between the two sets of results is caused by the different truncation errors of 
the respective perturbative expansions. 

In order to account for the deviations of the data from a perfect plateau, e.g., due to lattice 
artefacts and the truncation of the perturbative expansions we have then applied a more compli- 
cated procedure than a simple fit with a constant. In particular, one should take into account that 
residual lattice artefacts and truncation errors might conspire to produce a fake plateau. Therefore 
we consider it important to include correction terms for both types of errors. We have tried to 
incorporate higher terms in the perturbative expansions of AZ s (/j ip ) and ^,_ M Qw(/ip) treating 
the corresponding coefficients as additional fit parameters. Similarly we have attempted to correct 
for discretization effects by including a simple ansatz for lattice artefacts. Again, the parameters 
in this ansatz have to be fitted. Nevertheless, the number of fit parameters will not get too large 
because we fit the data for all four j3 values simultaneously. Only the quantities Z RGI (a), our final 
results, depend on j3, the other parameters do not. 

When we perform the fits we make the following choices. In the expansions for -^ri'.mom^p) 
and AZ 5 (/i p ) originating from continuum perturbation theory we use as many terms as are avail- 
able. The same applies to the /3 function used when computing the running coupling g s (/j, p ). We 

choose the MOM scheme as the intermediate scheme S and expand ^_ M qm * n ^ ne MOMgg 
coupling. All data for /x 2 > 10 GeV 2 are included in the fit. The correlations between the data 
at different momenta but the same f3 are not taken into account. Two examples of such fits are 
shown in Fig. [TT1 More details concerning the fit procedure can be found in Appendix [Ej 

In some cases, the data on our coarsest lattice (/3 = 5.20) are not very well reproduced by the 
fit, see, e.g., the upper plot in FigQTJ Excluding these data would, however, lead only to tiny 
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FIG. 11: Z KGl (perturbatively subtracted) for the operators Ot (upper plot) and V2 a (lower plot) as a 
function of the renormalization scale. Also shown are the fit curves used for the determination of Z . 

changes in the results. Therefore we have kept f3 = 5.20 in the fit for all operators. Note that the 
"divergence" of the fit curves (and the data) in the vicinity of fip = 3 GeV 2 is mainly caused by 

the Landau pole in the renormalized coupling constant in the MOMgg scheme. 

While the fits for the subtracted data are reasonable it was hardly possible to obtain a sat- 
isfactory fit for the unsubtracted numbers. Although plots of unsubtracted data do not differ 
dramatically from plots of subtracted data (see Fig. [12] for results for an operator where no sub- 
tracted data are available), the fit curves look quite strange. Therefore we have to conclude that 
our fit procedure is applicable only to perturbatively subtracted data and we must apply a differ- 
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ent procedure to unsubtracted data. So we choose the following method. We read off Z RGI at a 
reasonable value of /i 2 and take as the error the maximum of the differences with the results at 
one lower and one higher value of the scale. The choice of these three scales is to some extent 
dictated by the necessity to avoid large lattice artefacts as well as large truncation errors in the 
perturbative expansions. We take the values \? p = lOGeV 2 , 20GeV 2 and 30GeV 2 . 
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FIG. 12: Z KGl for the operator Oh 2 b as a function of the renormalization scale. 



XI. RESULTS 



Before we present and discuss our results we have to consider the influence of the two parameters 
that enter our analysis: the physical value of ro and roAj^g. The perturbative expressions which 
are needed in the evaluation of Z RG1 are functions of /ip/A 2 -^-, where ^ 2 is related to the momenta 



in lattice units q z by a 2 // 2 = q 2 . Since we use ro to set the scale we write 



(r /a)V 



(90) 



so that 



j4_ 

A 2 - 

MS 



(ro/a) 



2 „2 



(^Ams)' 



(91) 



This shows that the physical value of ro has an influence only on the scale \x v to which a particular Z 
value is associated. As the question of the scale at which perturbation theory becomes applicable is 
not completely immaterial we shall set ro to a reasonable value, for which we take 0.467 fm 48, 4^] 
(see also Ref. [io|). However, the precise number does not matter too much because ro enters only 
logarithmically. 
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The value chosen for roA^jg, on the other hand, has an impact on the results for the renormal- 
ization factors. In particular, varying roA^jg modifies the scale dependence of the right-hand side 
of Eq. ([541 an d improves or deteriorates the appearance of the plateau. We take roA^r = 0.617 
from Ref. 31], which is consistent with the value found in Ref. [511 ] . 

In order to estimate the systematic errors due to the uncertainties in the values of vq and roA^jg 
we have repeated our analysis using ro = 0.467 fm, roA^jg = 0.662 and ro = 0.5 fm, roAj^g = 0.617 
instead of our standard values ro = 0.467 fm, roAj^g = 0.617. Note that r Aj^g = 0.662 results 
from roAjyjg = 0.617(40) (21) as given in Ref. 3l| by adding the errors in quadrature. The ensuing 
differences will be shown in the same format as the (statistical) errors, however with the sign 
information included. For reasons of better readability they are given below the results themselves. 
The first (second) number corresponds to the difference caused by the variation of roA^jg (ro). For 
example, the entry 



0.45155(80) 
(568)(15) 

means that the analysis with ro = 0.467 fm and roA^jg = 0.617 produced the result 0.45155 ± 
0.00080 while using ro = 0.467 fm along with roA^r = 0.662 led to 0.45723 and working with 
ro = 0.5 fm, roA^jg = 0.617 gave 0.45170. Note that the first error ((80) in the example) is 
determined from the deviation of the Z RGI data from a perfect plateau as explained in more detail 
below, where also further sources of systematic errors will be discussed. 

In the cases where lattice artefacts have been subtracted perturbatively as explained in Sec. IIXI 
we have determined values for Z RGI by means of the fit procedure described in the previous section. 
The corresponding results will be called the fit results. 

For unsubtracted data we apply the alternative method mentioned at the end of Sec. [X] We 
evaluate Z RGl at the scales /i 2 = lOGeV 2 , 20GeV 2 and 30GeV 2 , interpolating linearly in // 2 
between adjacent data points. We take the value at 20GeV 2 as our central value and estimate 
the error from the maximum of the deviations of the values at the other two scales. This error is 
always larger (in most cases considerably larger) than the statistical error. These results will be 
called the interpolation results. 

Of course, the same method can also be applied to the subtracted data. Except for O at 
j3 = 5.20, the error of the interpolation results is again larger than the statistical error. So, for 
the operators for which perturbatively subtracted data exist we have interpolation results and fit 
results, both based on the subtracted numbers, as well as interpolation results extracted from the 
unsubtracted data. For the operators for which no subtracted data are available we have only the 
interpolation results. Note that the t>4 operators are particularly difficult: T-Q mri {p) vanishes on 
the diagonal of the Brillouin zone close to which all our momenta lie, and large lattice artefacts 
are obvious in the Monte Carlo data. So the corresponding nonperturbative results should be 
considered with caution. 

Finally, we have to decide which numbers we want to consider as the most reliable results to 
be used in the applications. It is clear that we make use of the perturbatively subtracted data 
whenever they are available. In these cases our fits seem to exploit the Monte Carlo data in an 
optimal manner, and therefore we take the fit results as our final numbers. However, the errors 
computed by the MINUIT program [5^] appear to be seriously underestimated as they are mainly 
determined by the statistical uncertainties. Hence we adopt a kind of hybrid approach taking the 
errors from the interpolation results (based on the subtracted numbers) because they take into 
account the deviation of our Z RGl data from a perfect plateau. The uncertainties due to the scale 
setting and the value of A^jg are again taken from the fit results. 

For the operators with two or more derivatives we do not have much choice. So we take the 
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TABLE V: Final nonperturbative results for operators with at most one derivative obtained with ro = 
0.467 fm and roAj^g = 0.617 from perturbatively subtracted data. Estimates of systematic errors have been 
included. 
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interpolation results (including the systematic uncertainties) as our final numbers. All our final 
results are collected in Tables [V] (operators with at most one derivative, based on perturbatively 
subtracted data) and [VI] (operators with more than one derivative, based on unsubtracted data). 

The impact of an uncertainty in ro or roA^jg is easy to quantify and therefore given in Tables IVl 
and IVII Further systematic errors are more difficult to control. In particular, the error caused 
by gauge fixing is hard to estimate reliably. However, as already remarked above, the existing 



investigations indicate that the "Gribov noise" does not exceed the present statistical errors [24 



26]. Since they are not the dominating uncertainty it seems justified to neglect the influence of 



Gribov copies although a more detailed study would clearly be desirable. 

The necessity of gauge fixing in the Rome-Southampton approach has also another consequence: 
The operators of interest, though of course gauge invariant, can mix with non gauge-invariant 
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TABLE VI: Final nonperturbative results for operators with two and three derivatives obtained with tq — 
0.467 fm and foA^jg = 0.617 from unsubtracted data. Estimates of systematic errors have been included. 
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(NGI) operators. However, in perturbation theory this effect shows up only in two-loop order and 
can thus reasonably be expected to be small. NGI improvement terms for the quark propagator 
were discussed at some length in Refs. 53j,|54|]. While 0(g 2 ) perturbation theory cannot distinguish 
between the gauge-invariant and NGI improvement terms, in an 0(g 3 ) calculation of the qqg vertex 
[55I ]. we could calculate the NGI improvement coefficient, and found it indeed to be numerically 
small. Mixing with gauge invariant operators, on the other hand, can in most cases be excluded 
by means of symmetry arguments (see Sec. HI)) . 

A few further systematic uncertainties can be estimated more easily. We have tested the sensitiv- 
ity to the chiral extrapolation by repeating the analysis employing a quadratic chiral extrapolation 
(see Eq. (l4"0")) ), This changed the results by less than 1%, except for the case of V4 at (3 = 5.20, 
where a change of 1.8% was observed. In order to estimate the error caused by the truncation of 
the perturbative series we have reduced the order of all perturbative expressions involved by one 
compared to the maximal value available. This led to changes of at most 1%. Finally, we have 
considered the uncertainty related to the chiral extrapolation of r^ja. Since little is known about 
the quark mass dependence of ro, we had to rely on some phenomenological ansatz [31] leading 
to chirally extrapolated values of r$/a with errors of the order of 1% (see Table ITTT)) . Varying the 
values of r$/a used in the analysis by 1% produced changes of at most 0.5% in the results for Z RGI . 
Thus it seems justified to assign an additional uncertainty of about 2% to our results. 

In the case of the perturbative estimates we consider the choices "bare PT", "TI PT", and 
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TABLE VII: Perturbative estimates for Z RGl based on one-loop lattice perturbation theory. The interme- 
diate scheme S is taken to be the MS scheme and roA^g = 0.617. In all cases Csw = 1 is used. 
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"TRB PT". Here "bare PT" and "TI PT" refer to the expressions AZ MS {fj, )Z^ c (fi ,a) pcvt and 
AZ MS (ii )Z^ e (fj, ,a)ti, respectively, both evaluated at /x = 1/a (see Eqs.£73]) and ([73])): "TRB 
PT" corresponds to the estimate by tadpole-improved, renormalization-group-improved boosted 
perturbation theory in Eq. ([82]) . Because of our choice /j,q = 1/a the perturbative estimates do 
not depend on tq. The bare PT and TI PT values do however depend on ^o^-ms as we ^ as 011 
the intermediate scheme S, which was taken to be the MS scheme in Eqs. (|73p and (|74p . We stick 
to this choice and set roA^jg- = 0.617. The TRB PT value, on the other hand, is independent of 
these choices. Note, however, that all three perturbative estimates depend on the chosen value of 
c sw . As remarked above, we have set c sw = 1. For the operators without derivatives we give the 
resulting numbers in Table IVIII The analogous results obtained from the two-loop calculations as 
described in Appendix [D] are given in Table IVIII1 All these perturbative numbers apply to any 
member of the H(4) multiplet to which the operator listed belongs. 

Let us now compare the results obtained by the various methods, i.e., by the different procedures 
of extracting Z RGI from the Monte Carlo data and by the different versions of lattice perturbation 
theory. In particular, for the operators for which perturbatively subtracted data exist we can 
compare the results extracted from the perturbatively subtracted data, both by interpolation and 
by means of the fit procedure, and the interpolation results based on the unsubtracted numbers. 
Of course, ideally they should agree within the errors. In reality, this is not always true. Note, 
however, that the errors of the fit results only account for the (rather small) statistical uncertainties 
of the raw data while the errors of the interpolation results are dominated by systematic effects. 

Figures [T3l 1141 and 1151 give an overview of our results for f3 = 5.40. The corresponding plots for 
the other (3 values look similar. For the operators without derivatives (see Fig. [T3]) the nonper- 
turbative results obtained with the different methods (with and without perturbative subtraction 
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TABLE VIII: Perturbative estimates for Z RGI based on two-loop lattice perturbation theory. The interme- 
diate scheme S is taken to be the MS scheme and r Aj^g = 0.617. In all cases the one-loop value for csw is 
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TI PT 


0.6943 


0.7013 


0.7068 


0.7204 




TRB PT 


0.7067 


0.7129 


0.7179 


0.7301 




bare PT 


0.8101 


0.8124 


0.8142 


0.8191 




TI PT 


0.7320 


0.7380 


0.7428 


0.7545 




TRB PT 


0.7401 


0.7456 


0.7500 


0.7608 




bare PT 


0.9789 


0.9857 


0.9905 


1.0009 




TI PT 


1.0007 


1.0064 


1.0105 


1.0190 




TRB PT 


1.0085 


1.0136 


1.0178 


1.0285 



of lattice artefacts) are well consistent in most cases. The one-loop perturbative estimates are 
larger, but tadpole improvement works. TRB perturbation theory, on the other hand, leads to 
further improvement only in a few cases, for some operators it is even worse than ordinary tadpole 
improved perturbation theory. 

In the case of the operators with one derivative (see Fig. fl~4"l) the agreement between the different 
methods for the nonperturbative results is less convincing, in particular, the interpolation results 
obtained without perturbative subtraction of lattice artefacts lie about 1 - 2% lower. Also the 
numbers from bare perturbation theory are smaller than the nonperturbative results. Again, 
tadpole improvement moves the perturbative estimates in the right direction, though too far in 
some cases. However, TRB perturbation theory leads to a significant overestimation. 

For operators with two derivatives (see Fig. [ToT) the only nonperturbative numbers we have at 
our disposal are those obtained without perturbative subtraction of lattice artefacts. Since the 
corresponding results for operators with one derivative lie consistently below those coming from 
the perturbatively subtracted data, it is tempting to guess that this is also the case for operators 
with two derivatives though we cannot quantify the difference. The behavior of the perturbative 
estimates is similar to that observed for operators with one derivative. Bare perturbation theory 
underestimates the results considerably, while tadpole improved perturbation theory comes much 
closer to the nonperturbative numbers. The results from TRB perturbation theory lie too high 
again. It should however be noted that the relative positions of the nonperturbative renormal- 
ization factors for almost all operators considered are surprisingly well reproduced by any of the 
perturbative estimates. 

In Fig. [16] we plot our fit results for the operators without derivatives together with the one- 
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FIG. 13: Results for operators without derivatives at /3 = 5.40. The filled symbols correspond to our fit 
results (circles) , interpolation results based on subtracted (squares) and unsubtracted (triangles) data. Our 
final numbers are the fit results with the errors taken from the interpolation results based on the subtracted 
data. The open symbols represent estimates from bare perturbation theory (circles), tadpole- improved 
perturbation theory (squares) and TRB perturbation theory (triangles) based on one-loop calculations. 

loop and two-loop perturbative estimates, again for (3 = 5.40. Let us comment on the numbers 
from bare lattice perturbation theory first, represented by circles in the figure. They exhibit the 
expected behavior: The two-loop results come closer to the nonperturbative numbers than the one- 
loop estimates, though only slightly in the case of the tensor current Oj^ v . Except for the tensor 
current, tadpole improvement works also in the two-loop approximation moving the perturbative 
values, indicated by squares, closer to the nonperturbative numbers. However, the results from 
TRB perturbation theory, shown by triangles, do not differ much from the values found by tadpole 
improved two-loop perturbation theory. At the moment it is unclear why the tensor current shows 
such a peculiar behavior. 

The perturbative estimates can easily be calculated at arbitrary values of the bare coupling 
constant g. However, for tadpole improvement one also needs nonperturbative values for uo (or for 
the average plaquette P = Ug) at these couplings. Such values can (approximately) be obtained 
from the results for P given in Table IIIII by a simple Pade fit taking into account the known two- 
loop expression for P [561 ] . In Fig. [T7] we plot the tadpole-improved perturbative results for the 
renormalization factor of the local axial current in the one- and two-loop aproximation along 
with our nonperturbative numbers and those from the ALPHA collaboration [571 ]. 

In a few cases we can compare our nonperturbative renormalization factors with results obtained 
by other methods. The renormalization factor of the local vector current OY, usually called Zy, 
can also be extracted from hadron three-point functions by considering the time component of the 
current and imposing charge conservation. Some time ago we have employed this approach in the 
case of the nucleon [58| on a subset of the gauge field ensembles used in the present work. The 
results are given in Table ILXl along with the renormalization factor of the local axial current Ofi, 
usually called Za, obtained by the ALPHA collaboration by means of the Schrodinger functional 
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FIG. 14: Results for operators with one derivative at j3 = 5.40. The filled symbols correspond to our fit 
results (circles) , interpolation results based on subtracted (squares) and unsubtracted (triangles) data. Our 
final numbers are the fit results with the errors taken from the interpolation results based on the subtracted 
data. The open symbols represent estimates from bare perturbation theory (circles), tadpole- improved 
perturbation theory (squares) and TRB perturbation theory (triangles) based on one-loop calculations. 
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FIG. 15: Results for operators with two derivatives at ft = 5.40. The filled triangles correspond to our 
nonperturbative results obtained by the interpolation method. The open symbols represent estimates from 
bare perturbation theory (circles), tadpole- improved perturbation theory (squares) and TRB perturbation 
theory (triangles) based on one-loop calculations. 
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FIG. 16: Results for operators without derivatives at /3 = 5.40. The crosses correspond to our nonper- 
turbative results obtained by fits of the subtracted data. The open symbols represent estimates from bare 
perturbation theory (circles), tadpole4mproved perturbation theory (squares) and TRB perturbation theory 
(triangles) in the one-loop approximation. The corresponding estimates based on two-loop calculations are 
shown by the filled symbols. 




FIG. 17: Renormalization factor of the local axial current as a function of g 2 . The curves represent one- and 
two-loop tadpole-improved perturbation theory. The circles are our nonperturbative results from Table |Vj 
the squares are numbers obtained by the ALPHA collaboration [57j . 
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TABLE IX: Results from alternative approaches. See Table 3 in Ref. [58j for Z v , Eq.(3J) in Ref. [59| for 



Z v (ALPHA), and Eq. (4.10) in Ref. |57j for Z A (ALPHA). The numbers for Z m in Ref. [60( are taken from 
Table 3 in that reference. 



Z RGI| 



=5.20 



7RGI 



^=5.25 



Z RGI| 



=5.29 



z RG1 \ B= 



/3=5.40 



Z V 

Z v (ALPHA) 
Z A (ALPHA) 
Z m [60j 
Z m (this work) 



0.7304(18) 

0.739(5) 

0.770(16) 

2.270(12) 

2.230(63) 



0.7357(13) 

0.744(5) 

0.774(13) 

2.191(24) 

2.222(62) 



0.7420(7) 

0.749(5) 

0.778(11) 

2.177(14) 

2.203(70) 



0.759(5) 
0.786(5) 
2.124(6) 
2.117(63) 



method [57f |. We also include the results for Zy from the ALPHA collaboration [591 ] . Note that in 
these determinations no gauge fixing is required. Hence the reasonable agreement of the numbers 
in Table HXl with those found in the present work indicates that the Gribov noise is small. 

In addition, we give in this table values for Z^ GI = (Za/Zp)^- gi , where Zp is the renormalization 
factor of the pseudoscalar density O p . The factor Z^ 1 renormalizes the quark mass as determined 
from the lattice axial Ward identity. In Ref. [(3(| we have calculated it by a method similar to that 
used in the present paper and applied the results in an evaluation of the strange quark mass. For 
easier comparison we also give the numbers following from Table [V] They differ from the older 
results by at most 2%. 

Factors for converting Z RGI to the MS scheme, i.e., AZ MS (^) _1 evaluated for our standard 
values of tq and r$ Atjc can be found in Appendix |F] 



XII. CONCLUSIONS 

As more and more detailed questions about hadron structure are treated in lattice QCD the 
renormalization of composite operators has become an important issue and perturbative as well 
as nonperturbative methods have been developed. In this paper we have presented results of a 
nonperturbative study in the RI-MOM scheme for a large variety of quark-antiquark operators, 
based on simulations with rif = 2 dynamical clover fermions. The results for the renormalization 
constants will be applied in the evaluation of phenomenologically relevant hadron matrix elements. 
Apart from these numbers, there are also a few lessons of a more general nature to be learned from 
our investigation. 

The renormalization factors connecting the bare operators on the lattice with their renormalized 
counterparts in some renormalization scheme, e.g., the MS scheme, depend on the cutoff used, 
the lattice spacing a in our case, and the renormalization scale M. The dependence on these 
two quantities should factorize, and this is indeed observed in a broad range of M. However, 
the available results from continuum perturbation theory for the anomalous dimensions and the 
(5 function can describe the M dependence only for relatively large values of the scale, above 
M 2 « 5GeV 2 . 

Only in the region where the scale dependence is well described by continuum perturbation the- 
ory is it possible to extract reliable values for the renormalization factors. On the other hand, for 
large values of the renormalization scale lattice artefacts may jeopardize the whole approach. It is 
therefore important to keep discretization effects under control, and we have seen that this purpose 
can be achieved (at least approximately) by subtracting lattice artefacts with the help of lattice 
perturbation theory. We did this at the one-loop level, but to all orders in a. Unfortunately, our 
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procedure turned out to be too complicated for operators with more than one covariant derivative. 
Alternatively, one can calculate the lattice artefacts proportional to a 2 in one-loop lattice pertur- 
bation theory. This has recently been done for operators without derivatives [61]. Since it should 
be possible to extend such calculations to more complicated operators, it would be interesting to 
see if subtraction of the a 2 contributions is already sufficient for our purposes. 

With the help of lattice perturbation theory one can not only calculate lattice artefacts, but 
also the renormalization factors themselves. However, due to the notoriously poor convergence 
properties of bare lattice perturbation theory some kind of improvement is mandatory, at least if 
only one-loop calculations are available. Indeed, we have found that tadpole improvement does 
quite a good job, although it is hard to predict how good the results actually are. For operators 
without derivatives, there are now even two-loop results. In bare perturbation theory they lead to 
a reduction of the difference with our nonperturbative renormalization factors, but the situation is 
less clear when tadpole improvement is included. Perhaps the ideal perturbative scheme is still to 
be found. 

Let us finally mention a few possible directions for future research. The RI-MOM scheme has 
the disadvantage that it requires gauge fixing. In principle, this problem could be overcome by 
working with correlation functions in coordinate space, and a first implementation of this idea has 
been published [62|. It seems, however, that very fine lattices are necessary for this method. 

Another possible modification concerns the choice of the momenta. In our application we have 
followed the original RI-MOM scheme, where the momentum transfer at the operator insertion 
vanishes. However, a generalization to nonexceptional momenta is possible 63h65|. 



A third variant of nonperturbative renormalization is motivated by the fact that the renor- 
malization condition (|28p involves only a particular trace of the vertex function T(p). On the 
other hand, we have the complete vertex functions (as 4 x 4-matrices after averaging over color) 
at our disposal, the bare ones computed nonperturbatively on the lattice as well as the renormal- 
ized ones calculated perturbatively in the MS scheme. So instead of introducing the intermediate 
RI' — MOM scheme by imposing (|28p one could directly compare the bare nonperturbative vertex 
function T(p) with the renormalized perturbative vertex function T MS (p) in the MS scheme. Up 
to lattice artefacts we should have 

T m (p) = (z^ Ic y 1 Z^ c T( P ), (92) 

where T MS (p) as well as the renormalization factors Z^ gxe and also depend on the renormal- 
ization scale \x. An analogous relation should hold for the quark propagator and Zr^ e . 

Of course, it is not to be expected that (|92|) is satisfied exactly: Not only lattice artefacts would 
spoil the identity, but also the truncation of the perturbative expansion. So one would have to 
develop some kind of fit procedure for extracting Z^ e from (|92]). In any case, it might be an 
interesting exercise to see how well (f92j) is fulfilled for our data. 
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Appendix A: Operator bases 

In this Appendix we list the operator bases we used when calculating the renormalization factors 
with the help of Eq. (151]) . 



For V2, a (representation t% 6 \ C = +1): 

O w , l< M <i/<4. (Al) 

For V2 } b (representation r\ , C = +1): 

Oll + 22 -C»33-044, 033-044, dl " ©22 • (A2) 

For r2 t a (representation t 4 , C = — 1): 

) 5 

{Hi/} 



OLa, 1<m<^<4. (A3) 



For r2,6 (representation 74 , C = — 1): 

5 u + Ol 2 -Ol 3 -Ol, Ol 3 -Ol, 0\ x -0\ 2 . (A4) 
For /ii ;a (representation r 2 , C = +1): 

2 <{^3} + <{W <{<W l<^i<^2<^3<4. (A5) 

(8) 



For (representation r{ , C = +1): 



^122 ~~ ^133 ) ^122 + C?33 ~~ 2 0144 , 

®211 ~~ ^233 > ^211 + ^233 ~~ 2 0244 , 

^311 "^322 > ^311 t ^322 ZL/ 344 > 

^411 — C*422 ; C4II + CJ22 ~~ ^Ojoo . 



(A6) 



For (representation rj , C 



C{122} - C{133} j C{122} + C{133} ~ 2 C{144} , 

^{211} - C>{233} > C{211} + C{233} _ 2 C{244} , 

C{311} - C{322} i C{311} + C{322} ~ 2 C{344} , 

0{411} - C{422} 7 C{411} + C{422} ~ 2 C{433} • 



(A7) 



For vsa (representation ri , C = — 1): 



For r3 (representation t 2 , C = +1): 

O 
O 



2 

C{234} , C{134} i C{124} i C{123} • (A8) 

(8) 



,5 

{122} 


~ ^{133} ' 


^{122} 


+ ^{133} 


~ 2 C{144} 


,5 

{211} 


— n 5 

^{233} ' 


^{211} 


+ ^{233} 


~ 2 C{244} 


,5 

{311} 


~~ ^{322} > 


^{311} 


+ ^{322} 


- 2 0{ 344 } 


,5 

{411} 


^{422} > 


^{411} 


+ ^{422} 


~ 20 {433} 



(A9) 



For 02 (representation C = +1): 

^{234} i ^{134} ' ^{124} ) ^{123} ■ 

For V4 (representation t[ 2 \ C = +1): 

C{1122} + C{3344} — C{1133} — C{2244} , 
C{1122} + ^{3344} + C{1133} + ^{2244} ~ 20{ 1144 } - 20{ 22 3 3 } . 

For h2, a (representation C = —1): 

^14{23} + ^24{13} + ^34{12} = ~~ ^^J{123} ' 
3Cl3{24} + ^14{23} + 3C , 23{14} + ^24{13} ~~ 2< ^U{12} = _9 C3{124} ~~ 3 C4{123} ' 
2 ^12{34} + ^13{24} + ^14{23} ~~ ^23{14} ~ ^24{13} = 3 ^'l{234} ~~ 3 ^2{134} ■ 

For /i2,6 (representation Tg 3 \ C = — 1): 

^l{122} ~~ ^1{133} + ^2{233} ' 
-2< ^1{122} ~~ ^1{133} + 3 ^'l{144} + ^2{233} ~ 3 ^2{244} ' 
_ ^1{133} + ®1{1U} ~~ ^2{233} + ®2{2U} ~~ 2 ^3{344} ■ 

For /i2, c (representation r^, C = 



>T 

13{32} 


+ ^23{31} 


~~ U 14{42} 


Ly 24{41} 


>T 

12{23} 


+ ^32{21} 


v - y 14{43} 


L/ 34{41} 


>T 

12{24} 


+ C , f2{21} 


Ly 13{34} 


^43{31} 


)T 

21{13} 


+ ^31{12} 


( - / 24{43} 


^34{42} 


|T 

21{14} 


+ ^41{12} 


^23(34} 


^43{32} 


>T 

31{14} 


+ ^IljlS} 


L/ 32{24} 


^42{23} 



For ^ 2 ,d (representation Tg 6 \ C = —1): 

^1211 ~~ ^1222 + Cl3{32} + C23{31} + ^14{42} + ^24{41} ' 
Cl311 ~~ ^1333 + Cl2{23} + ^32{21} + Cl4{43} + ^34{41} ' 



? 1311 ~~ ^1333 + ^12{23} + ^32{21} + ^ , 14{43} + ^34{41} ' 

Cl411 ~ Cl444 + Cl2{24} + ^42{21} + ^13{34} + C43{31} ' 

^2322 ~~ ^2333 + ^21{13} + ^31{12} + ^24{43} + ^34{42} ' 

^2422 _ ^2444 + ^21{U} + ^41{12} + ^23{34} + ^43{32} ' 



C3433 ~~ C3444 + ^31{14} + CJl{13} + ^32{24} + C42{23} ■ 

For the vector current (representation t[ 4 \ C = — 1): 

Ol, 1<^<4. 
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For the axial vector current (representation C = +1): 



O: 



1 < /i < 4. 



(A17) 



For the tensor current (representation t-[ , C = — 1 



O 



(IV 



1 < fj, < v < 4. 



(A18) 



Appendix B: Example for conversion factors 

We want to explain in detail, how the difference between Z^3_ MnM and ZHq M ar i ses > consid- 



ering the operator 



O 



(Bl) 



for / i/ as an instructive example. The operator V2 a belongs to this multiplet of operators. 

Using dimensional regularization, straightforward perturbation theory in 4 — e dimensions yields 
in the Landau gauge 



r(p) = \{i^Pu + iuV(i) 



;>> - + ln(4vr) - 7E - ln(//// 2 ; 
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(B2) 



up to terms which vanish for e — > 0. For QCD we have Cjr = 4/3. In this Appendix we restrict 
ourselves to one-loop order. Hence the coupling constant g can be identified with the bare cou- 
pling or with some renormalized coupling. Moreover, to this order and in the Landau gauge the 
quark wave function renormalization constant Z q is equal to one in all schemes of interest to us 
(RL - MOM, MS and MOM), see also Eq. (IC331) . So we can ignore it in the following. 
For the operator (jBip we have 



vPji) 



(B3) 



and the term proportional to p represents an additional structure which is not a multiple of the 
Born term. The renormalized vertex function in the MS scheme reads 



9 



T Ms(P) = r Born(p) + J^2 C F { ?Bom(p) 



Hp 2 /^) 
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and we get 



7 MS 
J dimrog 



1 + ife CF -lG + ln(47r) - 7E ! +0i(/l 



(B4) 



(B5) 



On the other hand, we can represent T(p) as a linear combination of rBorn(p) and PuPviif/p 2 )- 
Requiring that in the MOM scheme the coefficient of T-Q Qvn {p) be unity for p 2 = we find 



7 MOM 
J dimreg 



1 + 



167T 2 



C F 



-(^ + ln(4vr)- 7E ) + | 



0(g 4 



(B6) 



such that 



MOM 



(p) = r Bo rn(p) + 



^1^2 °F | r Bornb) " g hl(p 2 '/ ' fl 2 



\v^\ (B7) 
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Finally we obtain from Eq. ([29 



7 ri'-mom _ -, , 9_r< 

Z dimreg ~ 1 + Yo^2 F 

Now we can calculate 



8/2 , .„ , \ 31 4 p 2 ^ 2 , 
3U + ln(4 " ) " 7E J + ^ + 3^T^) 



+ 0(g 4 



(B8) 



7 MS ( 1 1\ 

^RI'-MOMvW 



r^MS 

dimreg 
^•RI'-MOM 
dimreg 



16vr 2 



9 3p 2 (p2+ p 2^ 



(B9) 



and 



7MS 2 / Qi ' 

'MS /,.\ _ ^dimreg _ -. . 5 ^ 61 \ , /->/„4a 



J MOM 



(A*) 



7MOM 

dimreg 



1 + 



(BIO) 



in agreement with the results given in Appendix O 



Appendix C: Results from continuum perturbation theory 

In this Appendix we collect the results from continuum perturbation theory that go into our 
computations. They all refer to n c = 3 colors and Landau gauge, but the number of flavors rif is 
left free. Note that we quote only the papers which give the results with the largest number of 
loops. 

We begin with the coefficients of the (5 function (see Eq. fl48j) ). In the MS scheme they are given 
by (see Ref. Q) 

A) = 11 -\n f , (CI) 
38 

ft = 102 -yn/, (C2) 

2857 5033 325 2 
ft = — -^ f + -nj, (C3) 

149753 nrnA . / 1078361 6508 \ 
ft - — + 35646 - (-162- + ) 



( 50065 6472 , \ „ 1093 , 
+ (^62- + ^r C3 ) n / + T29- n /- (C4) 



In the MOMgg scheme one finds 35] 

186747 1683, / 35473 65 , \ / 829 8 \ , 8 , 

B 2 = C3 C3 )n f C3 ni + -tit , (C5) 

H 64 4 Sd V 96 6 / V 54 9 / 9 

20783939 _ 1300563 _ 900075 _ / 2410799 _ 1323259 _ 908995 \ 

^ 3 ~ 128 32 32" ~ V 64 144 ^ 144 ^) Uf 

( 1464379 12058, 7540 \ „ /3164 64 \ „ 320 d 

+ (-648 ^ ~ ^ ) n f-{-W- ¥ C3 ) n f + -81 U f ' (C6) 

while /3o and ft are scheme independent in the Landau gauge. 

We now turn to the coefficients of the anomalous dimension in the MS scheme. Our conventions 
have been given in Sec. IVIIl see in particular Eqs. (|44p and (|45p . For notational simplicity the 
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superscript MS will be omitted. Note that we assume an anticommuting 75, so the coefficients 7$ 
correspond more precisely to the so-called naive dimensional regularization. 

While OY and 0$ have of course vanishing anomalous dimension, we find for O s and O p (see 
Refs. 0,0]) 



To 
71 

72 

73 



404 40 



fU32 , 320 
3" 



-2I!>S+|— + — Ca)n f + — n f 



V 27 
4603055 271360 



81 



+ 



10484 
243 



27 
1600 



9 



C 3 + 17600C5 + 

, 320 > A 2 
C3 + ^-C4 M + 



280 2 

/ 

/ 183446 68384 



(C7) 
(C8) 

(C9) 



V 27 
/66 
V243 



+ 



Cs - 1760C4 



36800 



C5 n f 



/664 128 , , , 
' Cs ) nj . 



27 



(C10) 



For 0^ u we have 69J 



70 
7i 
72 



3' 

724 104 

~9 27 
105110 1856 



-n 



81 



27 



Ca 



( 10480 320, . 

+ — Cs n f 



81 



9 



9 



(Cll) 
(C12) 
(C13) 



The operators V2 a , 0,; 2 6 , O r2 a and C r2 b have the same anomalous dimension. From Ref. [70(] we 
get: 



(C14) 
(C15) 
(C16) 





64 


7o = 


¥' 




23488 512 


7i = 


nt , 

243 81 1 




11028416 2560 


72 = 




6561 81 



/ 334400 2560 , . 

+ -r^Cs n f 



2187 



27 



1792 2 
729 / 



For O vs , O vz a , O a2 and O r3 we extract from Ref. [7l| : 





100 


7o = 


~9~' 




34450 830 


7i = 


243 81 n/ ' 




64486199 2200 


72 = 


26244 ' 81 



Cs 



/ 469910 4000 



2187 



+ 



27 



Cs n f 



2569 
729 



n f . 



The anomalous dimension of V4 can be found in Ref. [7 





628 




7o = 


45"' 






5241914 26542 




7i = 


30375 2025 


n f , 




245787905651 


11512 


72 = 


82012500 


405 



6- - 



/ 726591271 5024 



2733750 



+ 



27 



Cs n f 



384277 
91125 



-11 4 



(C17) 
(C18) 
(C19) 

(C20) 
(C21) 
(C22) 
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The three-loop anomalous dimension of transversity operators has been calculated by Gracey. For 



Oh x a and Oh x b we get from Ref. |72j 



70 = 8, (C23) 

71 = 124 -8n f , (C24) 

19162 /5608 320 \ 184 , _ . 

^ = —-{-2r + - c Tf-w n f (C25) 

In Ref. [z3] we find for O h , 2a , O h2b , O h2c and O h , 2d : 

70 = ^ , (C26) 
38044 904 



243 81 

17770162 1 280 / 552308 4160 \ 2408 2 

+ ^-ft-(w + ^r 6 J n/ "w"?' (C28) 



72 



6561 



Finally, we can take the anomalous dimension of the quark field to four loops from Ref. 74| : 



70 = 0, (C29) 

71 = — - g»V , (C30) 

20729 1100 40 9 , 

72 = - 79Cs - -g-n/ + -nj , (C31) 

2109389 565939, 2607 . 761525 / 324206 4582 320 \ 

73 = -gj _ 6 + _ a - _^- C 5 - (^-gj- + -^0. + 79Ci + n/ 

/7706 320 \ 2 280 , . 

+ "ST" + "S^s U + — nj . (C32) 



Let us now consider the coefficients needed for the conversion from the RI' — MOM scheme 
and the MOM scheme to the MS scheme, as defined in Eqs. (|56p and (I57p . respectively. Since the 
RI' - MOM scheme is in general not covariant, these coefficients may depend on the direction of 
the momentum p. In order to keep the paper at a reasonable length we refrain from giving the bases 
that are used in the representation of the vertex functions and enter the precise definition of the 
MOM scheme. The calculation makes use of the perturbative expressions for the vertex functions 
and of the ratio ^ M ™/^ are = ^ R bl~ e MOM /^ are . This ratio coincides with the quantity C 2 RI ' in 
Ref. [74] , where we can read off the expansion 

r^ lt j« ! + j« ! + .., (C33) 



with 



Z^S \ 16?H / \ 16tH 



■ ^ + 12(3 + ^/, (C34) 

439543 8009 . 79 ,. 1165. ( 24722 440 ,. \ 1570 
+ 



162 6 



I 79 1 165 / 24722 440 \ 1570 2 ,_ oe , 

-C3^ U Cs + C3 )nf ni. (C35) 

4 3 S \ 81 9 S J } 243 f V ' 
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Prom Ref. [74J we get for O s and O p 

ci = y, (C36) 
4291 152 83 

C2 = IT " T Cs " ¥ n/ ' (C37) 

3890527 224993 2960 / 241294 4720 80 \ 
° 3 ~ ^24 VT (3 + -9~ C5 " I ^43 27~ C3 + Y U " f 



+ 



With the help of Ref. [75j we find for 0^ 

ci = 0, (C39) 
3847 184 313 

c 2 = 1 C 3 -\ nt, (C40) 

1 54 9 s 81 n K J 

9858659 678473 1072 10040 / 286262 2096 80 \ 

C3 - "^9i^ + ^8^ Cs + ^r C4 " ^r Cs + {-^9 — 2^r C3 + t c v n/ 

/ 13754 32 \ 2 _ . 

- w + 8i <3 R (C41) 



For the vector and axial- vector currents 0^ and one finds from Ref. 75] 

ci = 0, (C42) 

C2 = f-^ + l"/) 12 . ( C43 ) 

ca = f-^f + 607C3 + - 32 C3 ) „, - f ^ * (C44) 



with i? = pft/p 2 . If the index /i is averaged over in the renormalization condition (see Eq. (|3ip ). R 
takes the value 1/4. 

Why don't C2 and C3 vanish although the vector current is conserved in the continuum? In 
the continuum the quark propagator S(p) and the vertex function T fJi (p) of the vector current are 
linked by the Ward identity 

ir» = ^-S-'ip). (C45) 

By Lorentz symmetry the massless inverse propagator must have the form 

S-\p)=\A(p 2 )f, (C46) 

where we expect A to depend logarithmically on p 2 / fi 2 . Therefore the vertex function has the form 

d^4 

r» = A(p 2 ) lfM + ^2prf , (C47) 
and the trace with the Born term gives 



dAv 2 c\A 
-tr ( 7 „I») = A(p 2 ) + 2p 2 — 2 p = ^(p 2 ) + 2p^J2 ^ ( ( MS) 
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The vector Ward identity therefore requires the existence of terms proportional to R in the trace, 
with a coefficient given by the logarithmic derivative of A. The coefficient c% vanishes because in 
the Landau gauge there is no term in the propagator of the form g 2 ln(p 2 //i 2 ), but in most other 
gauges ci is nonzero and we already have a term proportional to R at one loop. 
In the following cases we express the coefficients c\, c<i and C3 in the form 

Cl = c[ 1} + c (2) R, (C49) 
c 2 = 4 1} + 6 2 + 4 2) i?, (C50) 

C 3 = 4 1} + b 2 C ( l ) + 63 + (4^ + R > ( C51 ) 

where R contains the momentum dependence and is given in Table IXl The expressions R&(p) for 
j = 1,2,3,4 read: 



(n] = p 2 M-pI? pI{pI + pI-2p 2 , 

R (2) / s _ pIpIpI rv , 

n ^u\\r) 4. „2„2 , „2„2 V^' JJ ; 

[P) {pI+pI){pI-pD 2 + {pI + pI){pI-pD 2 [ j 



^ {4) (p) 



((p? + p1)(p1 + pD - 2 p\pI - 2 pIpIY 



p\(p\ +pI~ 2 pIY + pI(Pi +P\- 2p§) 2 + pUpI +P\- 2p!) 2 + P 2 4 (pI +pl- 2p 2 ) 2 

(C55) 



For V2a , V2 b , O r2a and O r2 b one extracts from Ref. 72] 



"2 
,(1) 



,(2) 
"1 

,(2) 



1 24 

^' (C56) 

68993 160 A 2101 ,^ r „. 

^ + ^ C3 + W^' (C57) 

451293899 1105768 > 8959 , 4955 , ( 8636998 224 640 \ 

1 C3 U C5 + C3 H C4 n f 

157464 2187 324 s 81 V !9683 81 27 V 1 

( 63602 256 \ 2 

-(w + 243 C3 J n /' (C58) 

J, (C59) 

2224 40 40 , n x 

~ir-T«'+T"" (C60) 



26244 243 



( 2 ) 136281133 376841 43700 > / 15184 1232 \ 9680 2 . 

■ - H (3 C5 + C3 )n f ni . (C61) 

' ' ■' Sd 81 s V 27 81 J 729 1 K J 
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TABLE X: Momentum dependent factors R. For the definitions of the lengthier expressions R^\p) see 
Eqs. (IC52|) - (IC551) . 



O O 



p 2 (pI+pI) 

{p\-(pI+pI+pDI^ 



°^'° r >-> 2p 2 (pi + (pl+P 2 2 +P 2 3 )/9) 

-9pl(pl-(pl+pl)/2y 



P 2 ((4p? + p\ + P 2 3 )pI + (Pi - (A + P§)/2) 2 ) 

-^p\p\p\ 

P 2 {{p\Pa) 2 + {P2PA) 2 + (P1P2) 2 ) 

m{p\-pI) 2 {pI-pI) 2 
p 2 ((p? +pI) (pI-pI) 2 + (pI+pI) (pI-pIY 
pjpj 

ip 2 ^ u pI+pI 



n n 

l - / f2,b 1 ^r 2 ,b n 2 



1 

6p 2 



{pI+pI-pI-pI? , (pI-p\? , {pI-pI? 



— — 9 

O O — 



(^Lb) + ^(p) + ^iLCp) + 4\Up) 



V3>a , O a2 -— (r^ 3 ( P ) + R (2 l(p) + + Rg } 4 (p) 

p(fl«(p) + i2W(p)) 



(9 



In the case of the operators a and Oh x b one obtains from Ref. [7 

4" = -y, (C62) 
2237 + | f3 + | n/) (C63) 



.(1) 



18 3 * 3 

1852993 97391 79 7060 / 306881 122 80 \ 

Cs - T C4 - ^C5 + ^C3 + y U ) n f 



+ 



432 108 
/1160 32 ^ 2 



; f> = 4 2 > = 4 2 ) =0 . (C65) 
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Similarly, Ref. |7J yields for O h2a , O h2b , O, 



l>2. 



and O h2 



,(2) 



218 
"~2T ' 

669202 452 4394 

^645~ + T^ Cs + W n/ ' 
1020141085 59050063 . 

+ 77777777 C3 



7679 



+ 



157464 
/ 98639141 
V 98415 



43740 " 324 
12712 1040 

C3 + 



1215 



27 



12434 
27 ^ 5 

/ 177970 416. . . 
" V"656T + 243 C3 ) 



.(2) 



.(2) 



0. 



The operators O vs , V3 a , O rz and O a2 require the coefficients 



214 

'~27~ ' 
4763093 

29160 



152 32363 

H Cs H «f • 

5 s 1944 ; 



8619089351 12125507 , 

+ 77777777: C3 



1574640 
/ 1227463 



10935 



8599 
324 ( 



2525 ^ / 1364405723 814 1000 A 

^T-CS + 777777777, + TT77/C3 + ^777-C4 



V 1574640 



135 



27 



400 

V 52488 + 243 



n.-, 



,(2) 



,(2) 



J2) 



4 

9' 

4432 



(C66) 
(C67) 

(C68) 
(C69) 

(C70) 
(C71) 

n f 
(C72) 

(C73) 
(C74) 
(C75) 

(C76) 
(C77) 

(C78) 

(C79) 
(C80) 

162000' v ' ^ C81 ' ) 

from Ref. (7^]. When using Gracey's results given in Refs. 0, l73l . l75j it is important to note that 
Gracey's RI' — MOM scheme is not the same as ours. Furthermore, Eqs. (A.l) and (A. 9) in Ref. [73| 
are not quite correct, but we hope that we have worked with properly rectified versions. Another 



56 . 50 

H Cs Tit , 

135 15 s 27 ; 

279011797 1717789 



131220 



2430 



9370 



( 1665047 28 \ 4210 2 

1 Cs nt H n f 

7290 5 s / ; 729 f 



extracted from Ref. [T3f] . For VA we find 



,(2) 



c (2) 

c 2 



J2) 



7214 
"675"' 

764724499 1756 

3645000 + 45 
282373048664443 



5655503 



243000 
796627067 



-TL{ 



39366000000 
/ 1160956742099 
h V 984150000 
1 

40 ' 



546750 



<3 



43507 



3208 , 1256 , 

+ -^777-C3 + ^777-C4 



135 



27 



1620 



Ci 



38398 

1167227687 
32805000 



+ 



2512 
1215 ( 



731129 



432000 
1047728166241 

9331200000 

51959 2 
:n f , 



23 119 
9CT 3 + 1200 Uf 



109467991 
+ 2916000 ^ 3 



13111 

~64T C5+ l 



( 232632277 755 . 

+ 777777C3 n, 



19440000 972 
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correction concerns Eq. (4.4) in Ref. [72], where the coefficient of TpNj in the a 2 contribution 
should be 468 and not 486. 

The corresponding coefficients c^, c' 2 and c' 3 in Eq. (j57[) for the conversion from MOM to MS 
are obtained by setting R = 0. So in the cases where the coefficients q are independent of R we 
have d i = c$ and the RI' — MOM scheme can be identified with the MOM scheme, at least to the 
order considered. 



Appendix D: Lattice perturbation theory to two loops 



In the two-loop approximation of bare lattice perturbation theory we have 

-2 / „2 \ 2 



7 MS 
J bare 



pert 



1 + 



9 



16vH 



-7 ln(a/i) + zi) + 



\16ir 2 



(hln 2 (afi) +h\n(a/j,) + z 2 ) ■ (Dl) 



For the currents the coefficients l±, l 2 , z\, Z2 can be read off from Refs. 0, 0] as functions of c sw . 
Note that z\ = — Cj?A in the notation of Sec. IVIIII 

Anticipating that we may want to expand in a coupling ^lat different from the bare lattice 
coupling <?, e.g., the boosted coupling g n (see Eq. ([71]) ). we express the MS coupling as a 
function of <?lat: 

— + 2^ In(oAi) - t^ AT + ( 2 7T A T¥ ln(a M ) - ^ sLvt + 0(s£at) • (D2) 



9- 



A1S 



5lat 



(16vr 



Here tf AT =ti~Pi = 1) 2), where the constants pi encode the relation between g and <?lat: 

1 1 



- 2 V\~V29 +0(g ). 
5lat 3 



(D3) 



For at = 5 one has pi = 0, hence i^ AT = ij, and the relation between and g takes the usual 
form|0-Si- 

For setting up tadpole improvement we need the expansion 



u 



l + n 



9 



+ r 2 



(-> 



2 \ 2 



16tt 2 ' ' z \l6n 2 

2 /„2 n 2 



+ 0(/ 



(D4) 
(D5) 



where = r 2 — 167r 2 rij>i. The coefficients ri and r2 can be found from Ref. [56|. For an operator 
with ri£) covariant derivatives the tadpole-improved two- loop expression for the renormalization 
factor then takes the form 



Z barc(Mti = 
.2 x 2 



1—nn 



1 + 



ffLAT 
16vr 2 



-70 ln(a^) + zi + (n D - l)r%) 



( 2 \ z 
J (/! In 2 (a/i) + (Z 2 + 16vr 2 p l7 o - (n D - l)r l7o ) ln(a/x) 

+ 4 AT + (n D - l)r£ AT + \{n D - l)(n D - 2)r\ + (n D - l)r x zA + 0(<?- 



(D6) 



LATJ 



where z\^- = z 2 — \§Tt 2 p\Z\. Of course, for the currents we have tid = 0. As in 
Sec. IVIIII the corresponding estimates for Z RGl are finally given by AZ MS (/j,o) Z^ e (fiQ, o)pert and 
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AZ MS (fj,o)Z^ c (fiQ, a)ti with /xq = 1/a. For the expansion parameter we take <7lat = 9n = s/uq. 
Therefore we have 

p ' = -^ = i c H- < D7) 

F 1287T 4 v ' 

In order to implement TRB perturbation theory at the two-loop level we start from Eq. (|78p . 
Using the three-loop expressions for 7 LAT and /3 LAT we obtain 

Z RGI = ( 2/ 3 fMr j ^ exp{F(^ AT /(16vr 2 ))} (D9) 
W = " " 

where 

A CiA = In 11 4- 

/?Q A) 



with 

ijLAT a LAT o flLAT., , o o LAT no oLAT LAT 

4 PoP 2 ^Po P 2 

= In (l + + ^|^x 2 ) (Dll) 



and 

1 / J^/ft - (/VA,) 2 ^ 
/ 2 (x) = _ arctan * — — . (D12) 

y/WZM/fo - (pi/A)) 2 V 2 + WM x J 

The explicit expression for 7^ AT has been given in Eq. (I80p . For the tadpole- improvement factor 
we make the ansatz 

ul~ n ° exp ( Cl /i( 5 £ AX /(167r 2 )) + c 2 f 2 (gt AT /(16n 2 ))) (D13) 
determining the coefficients c\ and c 2 such that 

exp (ci/i(^at/(16tt 2 )) + c 2 / 2 ( 5 2 AT /(16vr 2 ))) = u^' 1 + 0(^ AT ) . (D14) 
Then the final result in TRB perturbation theory reads 

2 - — 

Z^h = ul~ n - (2/3ofp) ^ exp {cx/^tA^ 2 )) + c 2 / 2 (^ AT /(16vr 2 ))} (D15) 



with 



LAT„, , o a „LAT no oLAT„,LAT 



C2 = 



2«LAT 



A ^lat , nfa. n 



Finally we have to decide how to deal with c sw . We insert the one- loop expression 

c sw = 1 + 0.2685885 2 + 0{g A ) (D18) 

in the above expansions and re-expand the result in the coupling constant. This means that we 
set c sw = 1 in the one-loop coefficients and the two-loop coefficients get additional contributions 
proportional to the one-loop coefficient in the expansion of c sw . 
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Appendix E: Fit details 



In this Appendix we give the details of the fits applied in the cases where perturbatively sub- 
tracted data are available. 

Working with the Moop approximation of the conversion factor Z^ 1 ,_ MOM (M) we write gener- 
ically 

4--momW = 1 + E y-^-) + R K M f S ' = ZKM) SS ' + Rl{M) ss ' (El) 

with the remainder R C ^{M) SS = 0(g s (M) 2i+2 ). The scheme S' chosen for the coupling in which 
^RI'-MOM * s expanded could be the MS scheme as in Eqs. (f56l) . ([57]) . ([58]) . but another option 

would be the MOMgg scheme. 

Unfortunately, in many cases the perturbative expansion of ^rj'_ MO m * s no * ver ^ we ^ conver " 
gent, in particular for S = S' = MS. However, if one chooses S = MOM, Z| turns out to be equal 
to 1 for O s , O p , QT u , O hl a , O hl b , Ha , O h2 b , O h2 _^O h2 d . Generally, the use of the MOMgg 
coupling seems to improve the convergence for S = MS as well as for S = MOM. Working with 
S = MOM instead of <S = MS seems to have the additional advantage that at least some of the 
effects of ^j/_ MO m are shifted to the factor AZ S via the anomalous dimension. In AZ S we can 
then exploit renormalization group improvement. 

Similarly, using the n-loop approximation of the f3 function and the anomalous dimension we 
express AZ S (M) as 



x' 



= Af(M)e R «W (E2) 

with Rf l {M) = 0(g s (M) 2n ). From the MS anomalous dimension one can compute the anomalous 
dimension in the scheme S to n loops, provided the conversion factor Zg (and the (3 function) is 
known to n - 1 loops, see Eqs. HMD - © for S = MOM. 

With the help of the above representations of ^'_ M om an< ^ AZ 5 we get from Eq. ([88]) 

Z!(^) SS ' + R C M SS ') ZZ; M ° M (^ a) = Z™\a)A s M~ l . (E3) 



In this relation as well as in Eq. (|88p lattice artefacts vanishing like a power of a have been neglected. 
For larger values of \x v this is not justified any more, even after perturbative subtraction of lattice 
artefacts. Therefore we write 

7 W'-MOM/„ n _ 7 RI'-MOM/ , \ A(rl 2, ,2\ fv > A \ 

Z bare a ) = Z bare i^P' a iMC ~ A{d /J, p ) (E4) 

subtracting the (remaining) lattice artefacts A from the Monte Carlo data Z^ r ~ MOM (/j,p, o)mc- Of 
course, A could be much more complicated than a simple function of a 2 fi 2 , but in the end we have 
to restrict ourselves to a polynomial in a 2 fi 2 anyway. So we use this simplified expression already 
here for notational convenience and obtain 

yet., yW-MOM Z RGI (q)Af (/z p ) 1 e ss , , o 2 \ f-p K \ 
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Given that Z^(^i p ) ss ' = 1 + 0(g s (/ip) 2 ) we approximate this relation by 
7 c (ll y>5' 7 Ri'-MOM r , ( n) _ Z^j^A^)- 1 e- RS M 

Z dHp) Z bare <W a JMC - 1 i nc(,. \SS' ^ A \ a Hp) ■ \^ b ) 



As we shall parameterize (R^) ss , and A in the following and fit the corresponding parameters, 
the above approximation should not be problematic. On the left-hand side we have our (possibly 
subtracted) Monte Carlo results for the renormalization factors, extrapolated to the chiral limit 
and converted to the intermediate scheme S using the Moop approximation. These numbers are 
fitted with the expression on the right-hand side, where the values of Z RGI (a) at our four values 
of a are the desired numbers. 



More precisely, we set 



K (Hp) = fl 9* (Ppr + h 9" W n+2 + • • • , (E7) 

RK»pf S ' = h /' (M 2£+2 + b2 9 S ' (Hp? e+4 + ■■■, (E8) 
A{a 2 n 2 p ) = gi a 2 f4 + g 2 {a 2 fi 2 p ) 2 + ■■■ , (E9) 

where b\, /1, ... are the (fit) parameters. Remember that in Sec. IIIII we have argued that 0(a) 
lattice artefacts are absent. In principle, the coefficients g±, g 2 , ■ ■ ■ could depend on the coupling, 
i.e., on /3. However, in the range of couplings we have at our disposal the variation of powers of a 
is much larger than the possible variation (logarithmic in a) of the coefficients. Therefore it seems 
justified to neglect this dependence and to treat gi, g 2 , ■ ■ ■ as constants. This might also help in 
disentangling lattice artefacts from truncation errors. 

For the operators considered here, the MS anomalous dimension is known to three loops, in 
some cases even to four loops, see Appendix [Cj Upon combination with the four-loop f3 function 
one can thus reach at least n = 3. The conversion factors Zr^-mom an< ^ ^MOM' 011 ^ ne other 
hand, are known to three loops in all cases. 

There are quite a few parameters that can be varied in the analysis: 



1. the intermediate scheme S, 

2. the scheme S' chosen for the coupling in which ^j,_ MOM is expanded, 

3. the orders of the perturbative expansions used in (jE6|) . i.e., the numbers n and £, 



4. the order Op of the perturbative expansion of the f3 function inserted in (|49h when computing 
the running coupling g s {M), 



5. the number of terms N\, N 2 , N a taken into account in the correction terms (|E7p . (|E8p . (jE9|) . 
respectively, 

6. the fit interval. 



Ideally, the results should be independent of all these choices. Moreover, one would expect that 
a significant deviation from the "continuum limit" curve (obtained by setting A = 0) appears only 
for Hp values where the data show lattice artefacts, e.g., in the form of a violation of the scaling 
property ([89]) . Unfortunately, the fits of unsubtracted data do not follow this expectation, and this 
is the main reason why we consider the corresponding results as unreliable and do not apply our 
fit procedure to these data. 

Our final choices are motivated by the following observations. The plateaus in Z RGI look better 
for the choice S = MOM than for S = MS. This may be due to the above mentioned fact that the 
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perturbative expansion of ^x'-mom seems to °e better behaved for S = MOM. For the scheme 
S' the choice S' = MOMgg seems to be favorable. For the operator a a comparison between 

S = S' = MS and S = MOM, S' = MOMgg is shown in Fig. (TDJ 

Not surprisingly, the maximal values for n and I lead to the best plateaus. For the number of 
terms taken into account in (|E7p . (|E8j) . and (|E9p only and 1 are reasonable choices. 

As already mentioned, the perturbative behavior, i.e., the plateau starts only at rather large 
values of the scale, typically around /^ 2 « 5 GeV 2 . So the lower limit /^in °f the fit interval should 
be at least 5 GeV 2 . However, the precise value does not seem to be too crucial, both fJ^ in = 5 GeV 2 
and = 10 GeV 2 look reasonable. 

So we arrive at the following choices. In the expansions originating from continuum perturbation 
theory which enter the fit formulae we use as many terms as are available, i.e., we take for n and I 
the largest values possible. The same applies to Op, so we set Op = 4. We choose S = MOM and 

S' = MOMgg. All data for // 2 > 10 GeV 2 are included in the fit. It is clear that the perturbative 
corrections (|E7j) and (lE8j) are hard to distinguish when inserted in (lE6j) and it does not make much 
sense to include both of them. So we set N\ = and N2 = 1. Furthermore we choose N a = 1. 
Hence we end up with six fit parameters: the four values Z RGI (a) along with the coefficients 61 
and g\. The data points are weighted by their statistical errors although the deviations from the 
fit curves are mostly of systematic origin. 

Appendix F: Going from RGI results to values in the MS scheme 

In this Appendix we collect the factors by which one has to multiply Z RGI in order to obtain the 
corresponding number in the MS scheme. They are given in Tables IXI) IXIII and IXIIII for various 
values of r$ and roA»-rq-. 



TABLE XI: Factors for converting Z KGl to the MS scheme obtained with r = 0.467 fm and roA-g^ = 0.617. 



Op. 


= 4 GeV 2 


^ 2 = 5 GeV 2 




1.40701 


1.44044 




0.91926 


0.91089 




1.04534 


1.04291 




u «2,.l ^V 2 ,bl ^2,«1 ^T2,b 


0.71544 


0.70183 


Oh!, a , O hlb 


0.69538 


0.68009 




0.58648 


0.56943 


Oh 2 , a , o h2:b , o h2c , o h2d 


0.57878 


0.56107 




0.50844 


0.49008 
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TABLE XII: Factors for converting 


; Z RGI to the MS scheme obtained with r 


= 0.467 fm and r A m = 0.662. 


On 


fi 2 = 4 GeV 2 


ii 2 — 5 GeV 2 

tx — U VJC V 




1.OO014 


1 A 1 G^O 

i.4iyoz 




Q94Q1 


u.y iouy 


A, 


1 .U^ ( U4 


i 0/1/1/11 


//~~\ rf\ 


0.72464 


0.71029 




0.70575 


0.68959 


^1> 3 J ^113,0! ^021 ^r 3 


0.59809 


0.58001 


fc 2 ,a> O h2b , O ll2c , O h2d 


0.59085 


0.57205 


a 4 


0.52100 


0.50145 



TABLE XIII: Factors for converting Z RG1 to the MS scheme obtained with r = 0.5 fm and roA^g = 0.617. 



Op. 


p? = 4 GeV 2 


p 2 = 5 GeV 2 




1.42764 


1.46022 




0.91406 


0.90608 


z q 


1.04382 


1.04155 


o o o o 


0.70698 


0.69402 




0.68587 


0.67134 


CO C9 

^d 3 j ^1)3,0' ^02 j 


0.57587 


0.55972 


^Va! ft 2 ,b) h.2,a) ^Vd 


0.56775 


0.55100 


a 4 


0.49699 


0.47969 
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